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ABSTRACT 



Observations indicate that massive stars in the Galaxy form in regions of very high surface density, 
£ ~ 1 g cm~ 2 . Clusters containing massive stars and globular clusters have a column density comparable 
to this. The total pressure in clouds of such a column density is P/k ~ 10 8 — 10 9 K cm -3 , far greater 
than that in the diffuse interstellar medium or the average in giant molecular clouds. Observations show 
that massive star- forming regions are supersonically turbulent, and we show that the molecular cores 
out of which individual massive stars form are as well. The protostellar accretion rate in such a core is 
approximately equal to the instantaneous mass of the star divided by the free-fall time of the gas that 
is accreting onto the star (Stahler, Shu, & Taam 1980). The star-formation time in this Turbulent Core 
model for massive star formation is several times the mean free-fall time of the core out of which the 
star forms, but is about equal to that of the region in which the core is embedded. The high densities in 
regions of massive star formation lead to typical time scales for the formation of a massive star of about 
10 5 yr. The corresponding accretion rate is high enough to overcome the radiation pressure due to the 
luminosity of the star. For the typical case we consider, in which the cores out of which the stars form 
have a density structure p oc r~ 15 , the protostellar accretion rate grows with time as cx t. We present 
a new calculation of the evolution of the radius of a protostar and determine the protostellar accretion 
luminosity. At the high accretion rates that are typical in regions of massive star formation, protostars 
join the main sequence at about 20M Q . We apply these results to predict the properties of protostars 
thought to be powering several observed hot molecular cores, including the Orion hot core and W3(H20). 
In the Appendixes, we discuss the pressure in molecular clouds and we argue that "logatropic" models 
for molecular clouds are incompatible with observation. 

Subject headings: hydrodynamics — ISM: clouds — stars: formation 



Massive stars are fundamental in the evolution of galaxies since they produce the heavy elements, energize the interstellar 
medium, and possibly regulate the rate of star formation. Remarkably little is known about how massive stars form, 
however: the problem is difficult observationally because massive star formation occurs in distant, highly obscured regions, 
and it is difficult theoretically because of the many processes that must be included. Even such a basic parameter as 
the time it takes to form a massive star has been uncertain. This time scale, or equivalently, the protostellar accretion 
rate, affects the luminosity of the protostar (particularly for masses to* < 1OM (Palla & Stahler 1992) and the strength 
of protostellar outflows (Richer et al. 2000). Arguments based on extrapolating from low- mass star formation lead to 
formation times t*t > 10 6 yr, a significant fraction of the main-sequence lifetime of the star (Bernasconi & Maeder 1996; 
McLaughlin & Pudritz 1997, hereafter MP97; Stahler et al. 2000). Comparison with observations of hot molecular cores 
(Osorio, Lizano & D Alessio 1999; Nakano et al. 2000) suggest substantially smaller time scales, i*/ < 10 5 yr. An analysis 
based on observations of protostellar outflows suggests i*/ ~ 3 x 10 5 yr (Behrend & Maeder 2001). The small (~ lx 10 6 yr) 
spread in ages of stars in the Orion Nebula Cluster (Palla & Stahler 1999), where there is no evidence that the higher 
mass stars have formed systematically later compared to the lower-mass population, sets an upper limit of t*f < 1 Myr in 
this case. What has been lacking is an adequate understanding of how the formation time is governed by the conditions 
in the gas out of which the star forms. 

Our understanding of low-mass star formation is on a far better footing, since it has received much more observational 
and theoretical attention (Shu, Adams & Lizano 1987). Low-mass stars form by accreting gas from a molecular "core" 
in which gravity overcomes thermal and nonthermal (magnetic and turbulent) pressure gradients. Shu (1977) considered 
the collapse of a singular isothermal sphere, finding 



where Cth is the isothermal sound speed and the numerical evaluation assumes tihc = 0.2nu 2 . Observed temperatures of 
10 to 20 K in regions of low-mass star formation imply accretion rates of about 10~ 6 to 10 -5 M Q yr^ 1 , consistent with 
the inferred values of £*/ for low-mass stars in these regions (Lada 1999). 

There are two difficulties in extending this theory to high-mass stars. The first, discussed in some detail by Stahler et 
al. (2000), is that the predicted accretion rate depends only on the temperature of the gas. Once massive stars form, 
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the gas may be heated to temperatures ~ 50 to 100 K, but the first massive stars that form in a region will emerge from 
gas at 10 — 20 K and will have low accretion rates and correspondingly long formation times. The second difficulty is 
feedback from the massive stars. Since the Kelvin-Hclmholtz contraction time is less than the accretion time for massive 
stars, they evolve along the main sequence while accreting. Massive protostars are thus very luminous, and it has been 
suggested that the radiation pressure and ionization they produce can halt the accretion and determine the upper limit 
of the stellar mass function (Larson & Starrfield 1971; Kahn 1974; Wolfire & Cassinelli 1987; Jijina & Adams 1996). This 
feedback is so strong that it is impossible to form stars as massive as those observed if the accretion is assumed to be 
spherical, and the discrepancy grows as the accretion rate is reduced. These considerations have motivated the radical 
suggestion that massive stars form via the coalescence of low-mass stars in order to achieve a more rapid build up of the 
final stellar mass (Bonncll, Bate & Zinnecker 1998). 

Recently, McKee & Tan (2002; hereafter MT) addressed the accretion-rate problem. The first step in resolving the 
problem of the apparently low accretion rates is to realize that, although equation (1) was derived for an isothermal gas, 
it should hold approximately when nonthermal support due to magnetic fields and turbulence is included as well (Stahlcr, 
Shu & Taam 1980; Shu et al. 1987). Observed cores have turbulent motions that increase systematically with radius 
(Larson 1981; Caselli & Myers 1995), and this leads to an increase in the accretion rate with time (Myers & Fuller 1992; 
Caselli & Myers 1995; MP97). Larger signal speeds allow for the hydrostatic support of denser gas cores, which then have 
shorter free-fall times and thus greater accretion rates once they become unstable. We term this the turbulent core model 
for massive star formation. 

The second step in resolving the accretion rate problem is the recognition that massive stars form in regions of very 
high pressure and density. MT showed that for typical pressures in regions of massive star formation (Plume et al. 1997), 
stars form in a time of order 10 5 yr. This result for the time scale is somewhat longer than that of Osorio et al. (1999), 
who inferred stellar masses and accretion rates by comparing calculated spectra with observations. 

The purpose of this paper is several fold. First, we determine the relation between the pressure in a molecular cloud 
and its surface density E. We show that observed regions of massive star formation, both Galactic and extragalactic, 
have E ~ 1 g cm~ 2 , corresponding to mean pressures P/k ~4x 10 s K cm~ 3 . Second, we extend the self-similar theory 
presented by MT to allow for magnetic fields and for a thermally supported core. Finally, we determine the radius and 
luminosity of accreting protostars, and use observed hot core luminosities to predict accretion rates and masses of several 
nearby massive protostars. 

2. THE PRESSURE AND SURFACE DENSITY OF REGIONS OF HIGH-MASS STAR FORMATION 

A cloud in hydrostatic equilibrium has a total pressure that is proportional to the square of the surface density, P ~ GE 2 , 
where £ = M(R)/ttR 2 (e.g., Elmcgreen 1989). The Appendix contains a detailed discussion of this relation, and for our 
fiducial case — including the effects of magnetic fields and allowing for embedded stars — we find that the mean pressure 
in a cloud is typically P ~ 0.88 GE 2 . What is the value of the surface density in regions of current and past massive star 
formation? 

Regions of high-mass star formation studied by Plume et al. (1997) are characterized by virial masses M v ; r ~ 3800 M Q 
and radii R ~ 0.5 pc. As discussed in the Appendix, the virial mass M v ; r is related to the actual mass M by the virial 
parameter, a v ; r = M v ; r /M . The virial parameter is about 1.3 — 1.4 for GMCs, whereas it is quite close to unity for cores 
that are actively forming stars. Although there is no direct determination of a v j r for massive-star forming clumps, Plume 
et al. (1997) regard the virial mass as the most accurate estimate for the mass (a v i r — 1)- Since this is consistent with the 
results for low-mass star formation discussed in the Appendix, we shall adopt a v ; r = 1 for our numerical estimates. The 
mean column density in the clumps studied by Plume et al. is then E c i ~ 1 gcm~ 2 = 4800M Q pc~ 2 ; the corresponding 
visual extinction is Ay — (Nn/2 x 10 21 cm~ 2 ) mag = 214E c i mag. From maps of a larger number of sources, Shirley et 
al (2002) find a median surface density (which they regard as a better characterization of their sample) from the virial 
theorem of E c i = 0.60 g cm~ 2 . Mueller et al. (2002a, b) have determined the gas masses of a number of these sources from 
observations of the 350 /xm dust emission. Using the Ossenkopf and Hcnning (1994) opacities, they find gas masses that 
are on average 3.4 times lower than the virial mass. The mean surface density E c i = 0.19 ± 0.12 g cm -2 for the sources 
they were able to model in detail, while their total sample gives a higher value, E c i = 0.73 ± 1.7 g cm~ 2 . Had they used 
the Pollack et al. (1994) opacities instead, the inferred masses and surface densities would have been larger by a factor 
of about 2.5. These data show that there is a considerable dispersion in the mean surface density of the sources. These 
column densities are far greater than those associated with GMCs (0.035 g cm -2 — Solomon et al. 1987) or regions of low 
mass star formation (the average column density in the C 18 cores in Taurus is 0.032 g cm~ 2 — Onishi et al. 1996). 

We now compare the column density observed in regions of massive star formation with the surface density of stars in 
both Galactic and extragalactic star clusters. Before doing so, we note that there are several effects that can alter the 
relation between the observed stellar surface density and that of the molecular gas out of which the cluster formed. First, 
since not all the gas in the cloud goes into stars, some of the gas is ejected, causing the cluster to expand. For example, 
if the star formation were 50% efficient and the gas were lost slowly, then the cluster would expand by approximately a 
factor of two from its initial size (Hills 1980). Including the mass lost from the cluster implies that the initial surface 
density of gas would have been 8 times greater than the final surface density of stars. On the other hand, if magnetic 
fields were important in the support of the cloud, then the stellar velocity dispersion would have been less than the virial 
velocity (Patel & Pudritz 1994), which would lead to an increase in the stellar surface density. If the size of the cluster is 
estimated from observations of the massive stars, then the resulting surface density will be larger than that of the total 
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star cluster insofar as the massive stars form preferentially in the inner regions of clusters (Bonnell & Davies 1998). This 
effect would be amplified if the cluster is old enough to have experienced significant mass segregation. Overall, we can 
expect the stellar surface density to be at best within a factor of a few of the initial gas surface density. 

First, consider Galactic star clusters. The virial mass of the Orion Nebula Cluster (ONC), including both stars and 
gas, is about 465OM (averaging the two models of Hillenbrand & Hartmann 1998). About half the mass of the cluster 
is within 0.8 pc, giving E ~ 0.24 g cm~ 2 . The Arches is a large cluster near the Galactic Center that contains many 
massive stars (Figer et al. 1999). Kim et al. (2000) have analyzed this cluster, taking into account the strong tidal 
effects experienced by a cluster so near the Galactic Center. They conclude that the initial mass of the cluster was about 
2 x 10 4 M Q . From their results, we estimate that the half-mass radius is about 0.4 pc, giving E ~ 4 g cm~ 2 . Figer et al. 
infer a very flat IMF for the Arches (dAC/dlnm* oc m~ o e vs. m~ 135 for the Salpeter IMF). If the actual IMF of the 
Arches is closer to the Salpeter value and the heavy extinction has led to an underestimate of the number of low-mass 
stars, then the surface density would be larger than this estimate. 

Globular clusters in the Galaxy no longer contain any massive stars, but they presumably did in their youth. There is 
a considerable dispersion in the properties of globular clusters; here we estimate the surface density of a typical one. The 
typical mass of a Galactic globular is about 2 x 10 5 M Q (inferred from Binney & Merrificld 1998). The data summarized 
by van den Bergh, Morbey, & Pazder (1991) give a median half-light radius for Galactic globulars of i? pro j.i/2 = 2.6 pc, 
where the subscript "proj" emphasizes that this is a projected radius. If the light traces the mass, then the spherical 
half-mass radius is R\/2 — 1.31i? pro j, 1/2 (Spitzer 1987), which is 3.4 pc. The surface density inside the spherical half-mass 
radius is then 0.5M/(7ri? 2 ^ 2 ) ~ 0.6 g cm~ 3 . The corresponding value of the pressure, ~3x 10 8 K cm- 3 , is similar to the 

estimates of Elmegreen & Efremov (1997) based on the same line of reasoning. 

It has been suggested that super-star clusters (SSCs) in external galaxies are globular clusters in the process of formation 
(e.g., Ho & Filippcnko 1996). The SSCs NGC 1569-A1 and A2 each have masses of about 4 x 1O 5 M (Gilbert & Graham 
2001) and projected half-light radii of about 1.7 pc (DeMarchi et al 1997), corresponding to E ~ 2.7 g cm~ 2 . A particularly 
dramatic SSC is the one observed by Turner, Beck, & Ho (2000) in NGC 5253: they estimate that a star cluster with 
total ionizing luminosity of ~ 4 x 10 52 s _1 is confined within a region about 1 x 2 pc in diameter. For zero age main 
sequence models (Schaerer & de Koter 1997), this corresponds to a total stellar mass of 0.6 - 1.5 x 10 6 M©, for Salpeter 
extrapolation to 1 and 0.1 M Q , respectively. If we take the projected half-mass radius to be 0.75 pc, then E ~ 20 — 50 g 
cm~ 2 . 

These results are collected in Table 1. It is striking that both the molecular regions where massive stars are forming 
and the star clusters that formed from them all have surface densities E ~ 1 g cm.- 2 (with the exception of the cluster in 
NGC 5253, which is somewhat larger). We first characterize the properties of molecular gas in such an environment and 
then determine the time scale for star formation. 

3. SELF-SIMILAR CORES AND CLUMPS 

Molecular clouds are inhomogeneous on a wide range of scales (Williams, Blitz, & McKee 2000), and numerical sim- 
ulations show that this is a natural outcome of supersonic turbulence, both with and without magnetic fields (Vazquez- 
Scmadcni et al. 2000). Following Williams et al., we define a star-forming clump as a massive region of molecular gas out 
of which a star cluster is forming; a core is a region of molecular gas that will form a single star (or multiple star system 



Table 1 

Characteristic Surface Densities of Regions of High-Mass Star Formation 



Object [Ref.] a 


M (M ) 


#1/2 (pc) 


E (gem- 2 ) 


P c l/fc (K cm" 3 ) 


Galactic Star-forming Clumps [1] 


3800 b - c 


0.5 C 


1.0 


4 x 10 s 


Orion Nebula Cluster [2] 


4600 b 


0.8 


0.24 


2 x 10 7 


Arches Cluster [3,4] 


2 x 10 4 


0.4 


4 


7 x 10 9 


Galactic Globular Clusters [5,6] 


2 x 10 5 b 


3.4 


0.8 


3 x 10 s 


NGC1569-A1,A2 [7,8] 


4 x 10 5 b 


2.2 


2.7 


3 x 10 9 


NGC5253 [9] 


0.6- 1.5 x 10 6 d 


1.0 


20-50 


2- 11 x 10 11 



References: (1) Plume et al. (1997); (2) Hillenbrand & Hartmann (1998); (3) Figer et al. (1999); (4) Kim et al. 
(2000); (5) Binney & Merrifield (1998); (6) van den Bergh et al. (1991); (7) Gilbert & Graham (2001); (8) DeMarchi et 
al. (1997); (9) Turner, Beck, & Ho (2000) 

b Virial mass estimates 

c The half-mass radius is not well-defined for the Plume et al. (1997) clouds, since the mass distribution on larger 
scales is not known. We therefore evaluate £ = M/ttR 2 using the typical radius and virial mass that they observe. 

"^Extrapolation from inferred LyC luminosity of H II region based on Salpeter IMF with a lower mass limit = 
1, 0.1 M Q . 
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such as a binary). Star- forming clumps are observed to be approximately gravitationally bound (Bertoldi & McKee 1992), 
whereas cores are necessarily gravitationally bound. 

Our basic assumption is that a star-forming clump and the cores embedded within it are each part of a self-similar, 
self-gravitating turbulent structure on all scales above that of the thermal Jeans mass (MT), an assumption that appears 
to be in accord with observation (Williams et al. 2000). Except on small scales, the turbulence is therefore necessarily 
supersonic, and we shall show that this is self-consistent below. We further assume that the structure is approximately 
stationary in time (i.e., it is not in a state of dynamical collapse or expansion); it follows that the star- forming clump and 
most of the cores within it arc in approximate hydrostatic equilibrium. However, some of the cores are gravitationally 
unstable and are collapsing on a dynamical time scale — these are the cores that are in the process of forming stars. Finally, 
we assume that the clump and the cores within it are approximately spherical. As shown by Bertoldi & McKee (1992), 
deviations from sphericity of a factor of a few are readily accounted for; in any case, observations of high-mass star forming 
clumps show that they are approximately spherical (Shirley et al. 2002). 

This model for regions of massive star formation is necessarily highly approximate. Regions of massive star formation 
are observed to have highly supersonic velocity dispersions (Plume et al. 1997) and are presumably turbulent. Large 
amplitude density fluctuations in such regions continually form, grow, and dissipate. Since we have assumed that a clump 
has a lifetime that is at least several times greater than its free-fall time, it follows that only a small fraction of the mass 
of the clump can be contained in density fluctuations that are gravitationally bound and are undergoing gravitational 
collapse — i.e., in cores that are in the process of star formation. The same conclusion applies to massive cores, which are 
themselves turbulent; if on the contrary, most of the mass of a core were in density fluctuations that were undergoing 
gravitational collapse, then the core would form a cluster of low-mass stars instead of a massive star. Calculations of 
turbulent, self-gravitating clouds indicate that indeed only a small fraction of the mass is gravitationally unstable at any 
given time (e.g., Vazquez-Semadeni, Ballesteros-Paredes, & Klessen 2002). 

In this section we shall work out the properties of this self-similar structure, leaving consideration of the effects of 
thermal pressure to §5. The results in this section apply equally well to an individual core or to a clump. For clumps, we 
allow for the possibility that only a fraction / g of the mass of the clump is gaseous, with the rest in stars. If stars are 
present, we assume that their spatial distribution is identical to that of the gas. We do not attempt to evaluate the mass 
distribution of the cores, which is related to the initial mass function (IMF) of the stars that form; instead, this mass 
distribution is assumed to be such that it leads to the observed IMF. Presumably this mass distribution is determined by 
the characteristics of the turbulence in the star-forming clump (e.g., Hcnrikscn & Turner 1984; Elmegreen & Falgarone 
1996). 

In a self-similar, spherical medium, the density and pressure each have a power-law dependence on radius, p cx r~ kp 
and P cx r~ kp . It follows that the sphere is a polytrope with P cx p 7p . In hydrostatic equilibrium we have (MP96; McKee 
& Holliman 1999) 

*p=9^-. k P = 7p k p = 22*-=2(k p -\). (2) 
Let c = (P/p) 1 / 2 be the effective sound speed. The equation of hydrostatic equilibrium then gives 

M=^, (3) 
where M is the total mass of the system, including any stars that are present. It follows that the gas density is 

P= ^GT 2 ' (4) 

where 

^^( 3 -M(fcp-i)/ g - 7p ^:ff /s (5) 

(cf. MP96). 

It is convenient to express the properties of the cores and clumps in terms of the pressure and mass. Since p — P/c 2 , 
equations (3) and (4) give 

_ { 2ttG 3 M 2 P \ 1/S 
C V Ak% ) ' 



/ AGM 



2nk 2 P P 



(6) 

1/4 



(7) 



and 

/ Ak 2 P P 3 \ 1/4 

In this subsection, A can refer either to clumps, which generally have / g < 1, or to cores, which are assumed to have 
/ g = 1; in the remainder of the paper, however, we shall use A to refer only to cores. 
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Observe that the only dependence of these properties on the polytropic index j p , or equivalently, on the power-law 
indices kp and k p , is in the factors A and hp, which are raised to relatively low powers. Studies of low-mass star formation 
often adopt the singular isothermal sphere as a model (Shu 1977), which appears to be approximately consistent with 
observation (Andre, Ward-Thompson, & Barsony 2000). Such a model has j p = 1, so that k p = kp = 2 and A = 1. 
No data are available on the structure of cores that are forming very massive stars. MP96 and MP97 have discussed an 
alternative model, the logatropc, in which the pressure varies as the logarithm of the density, but as discussed in Appendix 
B, this model is not physically realistic. Observations of molecular clouds show that the velocity dispersion a increases 
outwards (Larson 1981). Since we expect that <r oc c, and since c 2 = P/p oc r 2 ~ kp , this implies that k p < 2 and j p < 1 
(Maloney 1988). (If -f p is written as 7 P = 1 + 1/n, this corresponds to a negative index n; such a polytrope is therefore 
referred to as a "negative index polytrope.") 

For "high-mass" cores in Orion, Caselli & Myers (1995) find k p ~ 1.45 with a dispersion of ±0.2. According to van der 
Tak et al. (2000), the clumps in which high-mass cores are embedded have values of k p ranging from 1 to 2, centered 
around I.5. 1 Beuther et al. (2002) find k p ~ 1.6 ± 0.5 for their sample of massive star-forming regions, while Mueller et 
al. (2002b) find k p ~ 1.8 ± 0.4 for a sample of 31 sources (see Evans et al. 2002 for a review). For our numerical estimates 
of the properties of both high-mass cores and clumps, we shall adopt k p = 1.5, which corresponds to j p = 2/3, kp = 1, 
and A = 3/4. For comparison, the value k p = 1.8 found by Mueller et al. (2002b) corresponds to j p = 8/9, kp = 1.6, and 
A = 0.96. 

The value of 7 P determines the structure of a polytrope. As discussed by McKee & Holliman (1999), the stability of the 
polytrope depends on the adiabatic index 7 as well as on 7 p . If the polytrope is isentropic (7 = 7 P ), then the maximum 
density contrast between the center and edge of a stable molecular cloud is quite limited (< 14 for j p < 1). On the 
other hand, for non-isentropic clouds, an infinite density contrast is possible in a stable cloud provided 7 exceeds a critical 
value; for our fiducial case of j p — 2/3, this critical value of 7 is unity. The value of 7 for turbulent magnetic fields is 
4/3, significantly above this critical value. The effective value of 7 for the turbulent motions in molecular clouds is not 
known; this is a thorny problem, since the motions are inferred to damp very rapidly (Vazquez-Semadeni et al. 2000). 
The large density contrasts observed in molecular clouds suggest that, insofar as polytropic models are applicable, the 
effective value of 7 is large enough to render the clouds stable against dynamical collapse. 

3.1. The Mean Pressure in Clumps 

The mean pressure in a clump is P c \ = {l/V c \) J PdV. As discussed in Appendix A, this pressure is proportional to 
G£ 2 j, so we write 

P cl = foGZl, (9) 
where (j)p is a numerical factor of order unity. From equation (A5), we find 

<t>P = (j^J /g^geom^Savir, (10) 

where ^ gcom represents the effect of non-spherical geometry and 4>p = (c 2 )/(°' 2 ) the effect of magnetic fields. Here a is 
the one dimensional velocity dispersion; in the absence of magnetic fields, this is identical to the effective sound speed c, 
and 4>p = 1- Insofar as the clump can be approximated as part of a singular polytropic sphere, equations (A8) and (A9) 
show that the product tfipctvir is well determined. 

As discussed in the Appendix, clouds with aspect ratios of 2:1, either oblate or prolate, have gcO m — 1-1. Observed 
low-mass cores have aspect ratios of about 3:1 if they are oblate and 2:1 if they are prolate (Myers et al. 1991). Basu 
(2000) considered triaxial clouds and argued against larger aspect ratios for oblate clouds, which are formally allowed by 
the observations. Although the shape of high-mass cores is unknown, high-mass clumps appear to have quite spherical 
shapes (Shirley et al. 2002); we therefore adopt (f> gC om = 1 for our fiducial case. Observations suggest that the Alfven 
Mach number is of order unity, and in that case the discussion in the Appendix shows that <pp — 2.8. (By contrast, MT 
considered the non-magnetic case with <fip = 1-) The factor / g measures the fraction of the cloud mass in gas. Although 
the star formation efficiency averaged over all the molecular gas in the Galaxy is relatively small (~ 5% over the life 
of an association according to Williams & McKee 1997), we anticipate that the efficiency in the dense clumps forming 
massive stars could be considerably higher. We assume that a typical massive star- forming clump has about 1/3 of its 
mass in stars, corresponding to / g = 2/3. Finally, as discussed in the Appendix, the virial parameter a v i r is about unity 
in star- forming regions. Adopting these fiducial values for the parameters that enter <j>p, we find that the fiducial value 
is then <j>p ~ 0.88, so that P c i = 0.88GS 2 ! and 

Pd/fc = 4.25 x lO 8 !) 2 ! KcuT 3 . (11) 

1 van der Tak et al. (2000) remark that k p = 1.5 is just what is expected for regions undergoing supersonic infall. In our model, this region 
should not extend beyond the core radius — i.e., the radius of the region that is collapsing to form a massive star. The core radius has a fiducial 
value -Rcorc ~ 0.06 pc (see eq. 19), which corresponds to an angle of only a few arcscc at the typical distance of several kpc for observed 
massive star-forming regions. The beam sizes used to observe these regions are generally larger than a few arcsec, so the observed values of 
k p should typically refer to the clump in which the core is embedded rather than to the core itself. A central assumption of our model is that 
these star-forming clumps are in approximate hydrostatic equilibrium, which implies that the observed values of k p reflect the structure of this 
equilibrium, not a dynamical collapse. 
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3.2. Properties of Cores 

We are now in a position to express the properties of a core out of which a massive star forms in terms of the pressure 
of the clump in which it is embedded and the mass of the star that will form. The key point is that, on average, the 
pressure at the surface of a core will be the same as the ambient pressure in the clump there, Ps,corc = Pc\{r). We define 
the parameter </>p, C orc as the ratio of this pressure to the mean pressure in the clump. Since the clump is self-similar by 
assumption, we have 

-kp c 01. r^-^ r ^1 -2(*=p-l)/(3-fc P ) 



P c i(r) _ 3-fcp /_r ^ Kp _ 5 - 2k p 

' C " Pc.] " 3 Uc 



(12) 



For example, if the clump is a singular isothermal sphere (kp = k p = 2), then </>p. coro = \[M(r)/M c {\ 2 ; for our fiducial 
case of k p = 1.5, we have </>p, CO re = jW(r) /M cl ]~( 2 / 3 ). 

What is the typical value of (pp. COIC ? If the cores trace the mass of the clump, then for k p = 1.5 the value of <^p jC orc 
at the half-mass point (M = 0.5M c i) is 1.06. However, Bonncll and Davies (1998) have argued that, at least in the 
case of the Orion Nebula Cluster (ONC), the massive stars formed preferentially in the inner parts of the cluster. In 
particular, the Trapezium stars originated within about 30% of the half-mass radius of the ONC. If we assume that the 
ONC is similar to the Plume et al. clumps and adopt 0.3i?i/ 2 as the typical radius of formation of massive stars, then 
P Sj coro ~ P c i(0.3i?i/2)- For k p — 1.5, the pressure in the clump oc 1/r, so P SjCO re for massive star-forming cores is about 
3 times that at the clump surface (recall that we are identifying R\/2 with the radius of the Plume et al. clumps), which 
in turn is a factor (3 — fcp)/3 = 2/3 times the mean pressure inside the clump; hence, </>p C orc — 2 in this case. We shall 
adopt this as a fiducial value. 

The pressure at the surface of the core is then related to the surface density of the clump by 

P s , core 4*P, core <j> P GT^. (13) 

Numerically, we have 

P S ,coTe/k = 4.83 X lOVRcorc^pS 2 ! K cm" 3 -» 0.85 x 10 9 S 2 ! K cuT 3 , (14) 

where the final expression is for the fiducial case (<^p, C ore = 2.0 and <j>p = 0.88). 

To express the core properties in terms of the final mass of the star that will form out of it, m*/, we note that only 
a fraction e coro of the mass of the core will end up in the star; the rest is dispersed by the powerful protostcllar outflow 
associated with the star as it forms (e.g., Matzner & McKee 2000). We assume that the star formation efficiency 6 CO rc is 
constant throughout the accretion process, so that at any time the mass of a star, m*, is related to the mass of gas out 
of which it formed, M corc , by 

^* = ^core-^corc- (15) 

The value of e corc is uncertain for massive stars, since it depends on the unknown properties of protostcllar winds from 
massive stars. For low-mass stars, e coro <~ 0.25 — 0.75 (Matzner & McKee 2000). A similar analysis involving the scaling of 
low-mass outflows to high-mass protostars stars yields similar efficiencies (Tan & McKee 2002a, Paper II), so for numerical 
estimates we shall adopt e corc = 0.5 as a typical value. 

We begin with the one-dimensional velocity dispersion a = c/</>^ 2 . Since ^ is a number of order unity, a and c scale 
in the same way with r: 

a oc c oc r (2 - k ^/ 2 oc M< 2 - fc '>/ 3 < 3 - fc '>. (16) 
Equation (6) implies that the velocity dispersion at the surface of a core is given by 

/ , \ 1/8 / x 1/4 

gPcorcfflp \ ( m*f\ „l/4 



°-=w= lM \wt) \m) Kl kms • (17) 

-^fer^ tms -" <i8) 

where the last value is for the fiducial case. This velocity dispersion is considerably greater than the isothermal sound 
speed, c th = 0.27(T/20 K) 1 / 2 km s" 1 . As a result, the cores that form massive stars are necessarily supersonically 
turbulent and therefore clumpy. This is one of the principal results of this paper. It is often assumed that star formation 
can occur only when the gas becomes subsonic (e.g., Hcnriksen 1986; Vazquez-Semadeni et al. 2002), and this assumption 
leads directly to the problems with massive star formation that have been cogently summarized by Stahler et al. (2000). 
In our view, this assumption is unjustified, and we shall show that we can obtain a theory of massive star formation that 
appears to be consistent with the data by dropping it. 
The radius of a core is (eq. 7) 

A \ ( m* f \ 1/2 



i? corc = 0.050 2 — — ^- S i pc (19) 
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Finally, the density of H atoms at the surface of a core is (eq. 8) 

n H , s = 1.10 x 10 6 (AkULjl^pf 4 (i^) 1/2 ^f cm ~ 3 (21) 

i.n x 10 6 ( E 3 / 2 cm- 3 . (22) 

Note that equations for the clump velocity dispersion, radius, and density can be obtained from equations (17), (19), 
and (21) by making the replacements m*//e cor e M c i and A — > A corc / g ; by setting <fip lC ore = (3 — fep)/3 so as to set the 
pressure equal to its value at the surface of the clump (eq. 12); and by setting ow^s = c^sps so as to ensure that the 
clump is in hydrostatic equilibrium (eq. A8). 

How are the core properties related to those of the clump in which it is embedded? Consider first the surface density 
of a core, S C ore = ( m *// e corc)/^-Rcorc- Equations (7) and (13) give 

y /?P n 1/2 

^corc / ^i^p 



, Mcoro^pJ -» 1.22. (23) 

The cores thus have a column density comparable to that of the clump as a whole. Note that S c i includes the mass of 
the stars in the clump; the ratio of S core to the column density of gas in the clump is larger by a factor l// g . 

One of the important questions that must be addressed in any model for the formation of dense star clusters is how a 
large number of protostellar cores can fit into the small volume of the clump in which the stars are forming. This issue 
is automatically solved in our model since the cores are much denser than the clumps. We compare the mean density in 
a core, Pcoro to the mean gas density in the central part of a clump, p c \ = 4> p . corc p c i = 4> P , C orc3/ g M/(47ri? 3 ), where p c \ is 
the mean overall clump density and the factor <f) p . corc = (3 — k p )/3- (r/i? c i)- fep allows for central concentration. As before 
we assume the typical high-mass star-forming core is located at r = R c \/3 (Bonnell & Davies 1998), so that for k p = 1.5, 
4>p, core — 2.6. With the aid of equations (7) and (13) we find 



o , \ 3/4 , f M ^ 1/2 

3(PP, core \ ^ / ^core-^ 



Pel 4>p, core V 3 - kp J fV 4 \ m*f 



(24) 



,4i)"W'" (25) 



V4oooM y v m *f 

The mean density of stellar mass cores is about an order of magnitude larger than the mean gas density in the surrounding 
clump. Thus, even though the clump is forming a very dense star cluster, there is no "crowding problem." 

Finally, we compare the core radius with the tidal radius in the clump; for our model to be self-consistent, we require 
the core radius to be less than the tidal radius. If the clump were a rotationally supported disk, then the tidal radius is 
R t ~ 0.46[M coro /M(r)] 1 / 3 r, where r is the distance from the core to the center of the clump and M{r) is the clump mass 
within r (Paczynski 1971). For Galactic star-forming clumps, this is somewhat smaller than the core radius inferred from 

equation (7), i? C orc = fg 1 ^ 4 [M corc / M (r)] 1 / 2 r . However, molecular cloud cores in low-mass star-forming regions show little 
evidence for rotation (Goodman et al. 1992), and we assume that this remains true in regions of high-mass star formation. 
In that case, the clump is supported by turbulent motions and magnetic fields, and the tidal radius is determined by 



{k P 1) 



GM(r) 



Rt = (26) 



where we have assumed k p > 1. The ratio of the tidal radius to the core radius, 

r. fV 4 r M( r \ i 1 / 6 

(27) 



Rcoro {kp- 1)V3 [M cm 

is then somewhat greater than unity, as required. Note that this ratio increases weakly with the clump mass, so that the 
omission of effects associated with tidal gravitational fields becomes a better approximation for more massive star-forming 
clumps. 

4. SELF-SIMILAR HIGH-MASS STAR FORMATION 

Having determined the properties of massive star-forming clumps and the cores embedded within them, we now consider 
the rate at which a protostar grows once the core from which it is forming becomes gravitationally unstable. Dimensional 
arguments indicate that the accretion rate should be m* ~ M COTe /ts, where iff = (37r /32Gp) x / 2 is the free-fall time 
(Stahler, Shu, & Taam 1980; see also Shu et al 1987). Following MT, we write this as 

pi* 

m* — , (28) 
iff 

where 0» is a dimensionless constant of order unity. This equation has the same dependence on dimensional parameters as 
equation (1) if c t h is replaced by the virial velocity (Gm^/R) 1 ^ 2 (Stahler et al 1980; Shu et al. 1987), since (Gm*/i?) 3 / 2 /G oc 
m^(Gm^/R 3 ) 1 ^ 2 oc m*/ts- More precisely, equations (3), (4), and (28) imply 

4 c 3 

= 0*6 CO rc — ^kpA 1 ' 2 — , (29) 

7TV3 G 
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which is the generalization of equation (1) to the non-isothermal case. 

Equation (28) could be violated in the sense that </>* is much greater than unity only in the unlikely case that the star 
forms from a coherent spherical implosion; if the star formation is triggered by an external increase in pressure, </> t could 
be increased somewhat, but deviations from spherical symmetry in the triggering impulse and in the protostellar core will 
generally prevent <^>* from becoming too large. The equation could be violated in the opposite sense that 0* is much less 
than unity if the core is magnetically dominated, so that collapse could not begin until the magnetic field diffused out of 
the core. However, the magnetic energy in star-forming cores is generally less than or about equal to the gravitational 
energy (Crutcher 1999), so this effect is unlikely to be significant. Note that since the core out of which a massive star 
forms is highly turbulent, it will be clumpy and the accretion rate can be expected to show large fluctuations — i.e., in 
reality 0* will be a random function of time. 

MT showed that if the collapse is spherical and self-similar, then </>* is quite close to unity provided that the value 
of p entering tg is evaluated at the radius in the initial core that just encloses the gas that goes into the star when its 
mass is to*. The resulting value of 0* averages over the fluctuations that would occur in an actual turbulent core. For a 
polytropic sphere, p oc r~ fep and M oc r 3 ~ kp , which implies p oc M^ kp /^ kp ^ = M 2 /( 4 ~ 3 ' lp K Since we have assumed that 
the stellar mass is a constant fraction e coro of the gas mass (see eq. f 5), the accretion rate varies as 

to* oc m^p 1 ' 2 oc to*- 1/(4 ~ 37p) . (30) 

Integration yields 



t 



4-3 7p 



TO* = TO*/ . (31) 

where to*/ is the final stellar mass, which is attained at a time i*/ (MP97). This equation implies 

m* = (4-3 7p )^ = (4-3 7p )^(-M 3 3 \ (32) 
t U f V**// 

Note that for -f p < 1 the accretion rate accelerates (MP97); in particular, for the fiducial case 7 P — 2/3, we have to* oc t. 
From equations (28) and (32), the star- formation time, which is the time from the first formation of a protostellar core 
to the time the star reaches its final mass, is then 

(4_37p) 

t*f — iff, s , 

<p* 

where tg, s is the value of the free-fall time evaluated at the surface of the core — i.e., for the last parcel of gas to accrete 
onto the star. 

As discussed by MP97, termination of the accretion breaks the self-similarity once the expansion wave reaches to*/. 
This occurs at a time they label t cw , which is in the range (0.42 — 0.49)i*/ for < 7p < 1. If the boundary of the core is at 
constant pressure, a compression wave will be launched inward when the expansion wave reaches the boundary (Stahler 
et al. 1980). Equation (32) remains accurate until this compression wave reaches the center. Numerical calculations by 
Reid, Pudritz, & Wadsley (2002) show that in the logatropic case (which is related to the polytropic case with 7 P — > 0) 
this occurs at a time of about 0.7t*/. Thereafter, the accretion rate gradually decreases. 

MT used the results of MP97 to evaluate </>* in the non-magnetic case, which we label 0* non : 

0*non = (4 - 37 P )^, (34a) 
U f 

"(2-7 p ) 2 (4-3 7p )( 7 - 6 »/ 2 to ^ 1/(4 - 3 ^ ) 



7^3 



(34b) 



8 (5-3 7p )/2 

(see Fig. I), where to (7 p ) is a parameter that is evaluated by MP97. 2 For example, for the singular isothermal sphere, 
we have to = 0.975 and 0* non = 0.9757rV3/8 = 0.663. For other values of 7p in the range < 7p < 1, which is equivalent 
to I < k p < 2, the approximation 

</>* non ~ 1.62 - 0A8k p (35) 

is accurate to within about 1% (MT). For the fiducial case (k p = 1.5), we have ^ tnon = 0.90. We conclude that to* ~ m*/ts 
to within a factor 1.5 for non-magnetic, spherical cores in which the effective sound speed increases outward. 

We can estimate the effect of magnetic fields on the accretion rate from the work of Li & Shu (1997), who considered 
collapse of self-similar, isothermal, magnetized, toroidal clouds. The equilibrium surface density is £ = (l+Ho)c 2 h /(TrGzu), 
where vj is the cylindrical radius and Hq is a parameter that increases from zero as the magnetic field is increased. They 
show that the accretion rate is to* = 1.05(1 + Ho)c 3 il /G, which is larger than that in equation (1) by about a factor 
(1 + H ). However, equation (28) predicts to* oc Mp 1 / 2 oc M 3 ^ 2 /ru 3 / 2 oc Y?I 2 -U3 Z I 2 oc (1 + H f/ 2 . To reconcile this result 
with the correct value, we require 

+- * JTThoY^' (36) 

2 MP97 assumed that the gas was isentropic — i.e., that the adiabatic index 7 is the same as the polytropic index ~f p . As discussed above, 
the large observed density contrasts in molecular clouds suggest that 7 > 7 P . No calculations are available for this case, however, so we have 
adopted the isentropic results of MP97. 
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Li & Shu regard Hq ~ 1 as typical for low-mass star formation, which would lead to a reduction in the accretion rate by 
a factor of about 1.4. 

The star formation time t*f increases from 1.5(1 + HoY^ 2 tft,s to 3.5(1 + HoY^ 2 tfi,s as 7 P decreases from 1 to 0; for 
the fiducial case, t*f = 2.22(1 + HoY^ 2 ts.s- (Note, however, that models with j p close to zero are not realistic because 
they have no pressure gradient: kp — > from eq. 2. Also, see the discussion of logatropes in Appendix B). The 
ratio of the star formation time to the mean free-fall time of the core, iff, core = (37r/32Gp coro ) 1 / 2 , is larger by a factor 
(Pcorc/Ps) 1 / 2 = [3/(3 — /Cp)] 1 / 2 . For example, for k p = 1.5, we have t*f — 3.14(1 + H Y^ 2 tg^ coro . Finally, we note that if 
we use equation (25) to normalize the star formation time to the mean free-fall time of the clump in which the core is 
embedded, we find 

Hence, for typical conditions in Galactic regions of high-mass star formation, the time to form a star is about equal to the 
mean free-fall time of the parent clump, t*f ~ iff, c i- In super star clusters, which have M c \ ~ 10 6 M , the star formation 
time is somewhat less than this. 

4.1. Evaluation of the Protostellar Accretion Rate and Star Formation Timescale 

The general expression for the accretion rate in terms of the surface pressure and the final stellar mass can now be 
inferred from equations (31) — (33): 

m*=^(^Y, (38) 

where 

In the fiducial case (k p — 1.5), we have j = 0.5; for k p = 1.8, we have j = 0.25. Using equation (8) to evaluate 
ts,s « (1//Ocorc, s ) 1/2 , we find 

*• " * »*M*«i*& (£fc) m (£f£*f fe)' «o y t -, (40) 

i^T^i^T'^'- • (4i > 

where again we have used the symbol — > to emphasize that the final expression applies only in the fiducial case. Although 
the exponent j in the relation oc mi changes with k p , the numerical coefficient is constant to within 10% for 1.3 < 
k p < 2. The evolution of the mass-accretion rate and the mass is illustrated in Figure 2. Observe that our expression 



Table 2 

Model parameters and their fiducial values 



Parameter 


Definition 


Fiducial Value 


A 


(3-fc p )(fc„-l)/ fl 


1/2 (clumps), 3/4 (cores 


face 


fraction of core's gravitational energy radiated by star+disk 


0.5 


h 


gas mass fraction: M g /M — S g /S 


2/3 (clumps), 1 (cores) 


Ho 


magnetic field parameter in Li-Shu cores 


1.0 


kp 


clump and core density structure: p oc r~ kp 


3/2 


k P 


clump and core pressure structure: P oc r~ kp 


2(fc p - 1) = 1 


niA 


Alfvcn Mach number 


1.0 


Glvir 


virial parameter, 5(a 2 }R/(GM) 


1.0 


^core 


core star formation efficiency 


0.5 


7 P 


polytropic index: P oc p lp 


kp/k p = 2/3 


0geom 


deviation from sphericity: Pi^/ (R 2 Z) 
(c 2 )/(a 2 ) = 1.3 + 3/(2m 2 4 ) 


1.0 


(t>B 


2.8 


4>p 


P c i/(GS 2 i) = (3V2O)/<A com 0ija v i r 


0.88 


core 


Pci(r)/Pci 


2.0 


(fop, core 


Pc\(r)/p c \ 


2.6 


^*non 


m^ts/ra^ (non-magnetic case) 


0.90 


</>* 


0*„o„/(l + i?o) 1/2 


0.64 
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1 1.2 1.4 1.6 1.8 2 



Fig. 1. — Variation of model parameters (mo, </>*non) and results (m. t f, t,f) with k p . Values of mo are from Table 4 of MP97. Over the 
entire range of -y p and k p relevant to molecular clouds (0 < 7 P < 1, 1 < k p < 2), </>* n on — 1-62 — 0.48fc p to within about 1%. The star 
formation time decreases from 3.5% to 1.5iff as 7 P varies from to 1. The variation of m t f and t„j relative to the k p = 1.5 case is also 
shown. Note that the singular polytropic model in hydrostatic equilibrium breaks down for k p = 1, 7 P = 0, since then the pressure gradient 
vanishes (fcp =0). 
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for the accretion rate depends on two dimensional parameters, the final stellar mass, m*f, and the surface density of the 
clump, £ c i (or, equivalently, the pressure in the clump), as well as on the ratio of the current stellar mass, m*, to the final. 
We have expressed our results in this form in order to facilitate comparison with observation. In fact, the accretion rate 
depends on only two parameters, the enclosed mass, M, and the parameter K p = P/p 1 " , which is constant for polytropes 
(Tan & McKee in preparation). 

The fiducial case includes the effect of a magnetic field. How does the field affect the accretion rate? From equations 
(10), (13), and (38) we have 

,3/8 

ro, cx 0* P s 3 / c 8 orc cx 0* 4 /8 cx + ^ o)1/2 - 1-04, (42) 

since 4>b = 2.8 and Ho = 1 in the fiducial case. There is therefore very little difference between this case and the field-free 
case considered by MT. It should be kept in mind that we have assumed that the value of the virial parameter a v i r has 
been taken from observation and is therefore fixed as the field is varied. As discussed in Appendix A, if the value of the 
virial parameter is taken from theory, then a v ; r 0B is fixed and the accretion rate would vary as (1 + Ho)^ 1 ^ 2 - 
The corresponding value of the star-formation time is 



W = 3.15x10' ■ - > -^fH ' ' ,. V, (43) 




- 1M ' 10 '{mWj E <. (44 > 

For example, with the fiducial values of our parameters, a 100 Mq star forming in a clump with E c i = 1 g cm -2 has a final 
accretion rate of 1.1 x 10~ 3 Mq yr -1 and a star-formation time of 1.74 x 10 5 yr. If the value E c i ~ 0.24 g cm~ 2 inferred 
for the Orion Nebula Cluster was appropriate for the molecular gas out of which the Trapezium stars formed, then the 
time scale to form a 30 Mq star in Orion was about 3.8 x 10 5 yr. The dependence of t*/ on k p is shown in Figure 1. For 
example, had we adopted the value of k p — 1.8 found by Mueller et al. (2002b), the coefficient in equation (44) would 
have dropped to 8.9 x 10 4 yr. The star-formation time as a function of final stellar mass is shown in Figure 3. Effects of 
changing the pressure and using a different, though still representative, value of k p are also shown. 

Our results for the accretion rate and star formation time have been presented in terms of the surface density of the 
star-forming clump. However, if higher resolution observations are available, then it is convenient to express these results 
in terms of the measured density profile, bh = n-n yIC {(r/r IC {)~ kp . Note that the resolution must be high enough so that 
the data apply to the core forming an individual star, and not to the clump in which many stars are forming; for the 
fiducial case, this requires measuring the structure on a scale less than about 0.05 pc~ 10 4 AU (eq. 19). Expressing the 
density in terms of the mass, we have 

3/fe p c, n*=p/(3-fep) 
4 ^H^H.rcf r rof 



(3 - k p )M 



(45) 



where /xh = 2.34 x 10 24 g is the mean mass per H. We now ask, how long will it take the gas in M to form a star of 
mass e COIC M? Equation (33) yields 

■(a-*,) 6 "*' ' ^ n V2(3 " M 



1 / 3tt 



b*k p V8G 



1/2 



("H.rofMH) 3 V 47Te CO rer 3 ef 



(46) 



4 /10 8 cm- 3 \ /1000AU\ 3/2 / m* \ 1/2 , N 

-» 5.79 x 10 4 — yr. 47 

Had we chosen k p = 1.8, the typical value found by Mueller et al. (2002b), the result is only slightly different: 

on^ n4 /10 8 cm- 3 \ 5/4 /i000AU\ 9/4 / m* \ 3/4 
^ f ^8.07xl0 4 hmri- y r - 48) 

The median value of nn.ref a t r rof = 10 3 AU in the Mueller et al. sample is 2.4 x 10 7 cm~ 3 (based on the Ossenkopf & 
Hcnning 1994 opacities), which gives t,/ ~ 5 x 10 5 yr for a 30 M© star. Recall that this choice of dust opacity led to 
masses that averaged 3.4 times less than the virial mass; if the clouds were at their virial masses, the densities would be 
3.4 times higher and the median time to form a 30 M star would be 10 5 yr. 

Recently, Doty et al. (2002) have estimated the timescale for massive star formation using chemical clocks. From an 
analysis of the source AFGL 2591, they infer that the chemical timescale is in the range 7 x 10 3 — 5 x 10 4 yr, with a 
preference for 3 x 10 4 yr. The protostellar models discussed in §6 below suggest that the central star has a mass of about 
13 Mq. The time scale measured by Doty et al. refers only to a part of the star formation time, since the chemistry is 
strongly affected by the rapidly rising luminosity. The time scale they infer is in qualitative agreement with our results. 
Unfortunately, our results do not apply to the the case k p = 1 they used to model their data, since such a model cannot 
be in a self-gravitating equilibrium (see Appendix B); as a result, we cannot make a quantitative comparison. 
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Henriksen (private communication) has pointed out that the rn 1 /^ 4 dependence of the star-formation timescale was 
presented in Henriksen & Turner (1984) and Henriksen (1986). This scaling is to be expected, since it follows directly 
from the assumption of hydrostatic equilibrium. In these models it was assumed that stars form when the turbulence 
becomes sub-sonic and that the pressure is a constant in GMCs. As we have seen, the former assumption is incorrect for 
massive stars, and the latter assumption leads to a significant over-estimate of the star-formation time for massive stars. 

Adams (private communication) has pointed out that results similar to ours are implicit in Adams and Fatuzzo (1996). 
By dividing their estimate for the mass of a protostar by the accretion rate c 3 /G, one finds that the star formation times 
are about 10 5 yr and depend only weakly on the mass. Their work is based on the assumptions that (1) protostellar 
accretion is halted when the rate at which mass falls directly onto the star, rather than through a disk, is the same as the 
mass loss rate in the protostellar outflow; and (2) that the wind luminosity is a fixed fraction of the stellar luminosity. 
The validity of these assumptions, particularly for massive star formation, remains to be established. 

4.2. Turbulent Core Accretion vs. Competitive Accretion 

An alternative model for massive star formation is the competitive accretion model, in which massive stars form by 
accreting gas that was initially not gravitationally bound to the star. The accreting stars compete for the mass available 
in the clump, and it has been suggested that this process leads to the IMF (Bonncll et al. 1997; Bonncll et al. 2001). By 
contrast, the turbulent core model is based on the premise that most of the gas that goes into a star was initially bound. 
Since the core is turbulent and is embedded in a turbulent medium, this is undoubtedly an oversimplification — at the very 
least, some of the gas that was bound to the core at the onset of star formation will become unbound during the course of 
the star formation, and similarly some gas that was initially unbound will become bound. First we estimate the amount 
of mass exchange between the core and the clump, and then consider the rate of accretion directly onto a star once it has 
formed from its core. 

In general, accretion occurs at a rate given by 

m acc = ir p c \v rc \R 2 acc , (49) 

where p c \ = (p Pi C orcPci is the density of the ambient clump material, v rc \ is the relative velocity between this material and 
the accreting object, be it core or star, and i? acc is the accretion radius, which is the largest possible impact parameter of 
approaching gas that ends up bound to the core or star. Recall that the factor p . coro ~ 2.6 (§3.2) allows for the central 
concentration of clump density distribution and massive star formation. We assume that the relative velocity of the core 
or star and the ambient medium is ^/3cr c i = V3((/>p iCor c/</y corc) 17 ' 2 ^ 2 !). The rms velocity dispersion in the clump, (cr^), 
is related to the virial parameter a v - u by equation (A3), under the assumption that the clump is spherical. We then find 

- - * (fe) " = " 7 fe) "' El,v " (*Sk) " 4 km (50) 

4.2.1. Accretion by the core 

Equation (37) shows that in typical Galactic star-forming clumps, the star formation time is comparable to the free-fall 
time of the clump, which in turn is comparable to the time for the core to cross the clump. Since the surface density of 
the core and the clump are similar (see eq. 23) we anticipate that the core will interact with a mass comparable to its 
own as the star forms, or equivalently the rate at which a core interacts with clump mass should be similar to the rate of 
star formation. 

To make this idea more quantitative, we estimate i? acc for the core by considering particle orbits that intercept its 
surface. In this case i? acc = (1 + v^/v^Y^Rcore, where v csc is the escape velocity from the core, 

«~ = ™ ( " 4 ^ ( ^"—M " ! Z'J' km „-', (51) 



1/4 



The gravitational enhancement of the geometric cross-section is thus 

1/2 / Mr 1 : 



i + O^a-l-^f ] 1—^1 (53) 



SOMqJ V4000M Q 

The enhancement is significant only for massive cores in relatively low mass clumps. The rate of accretion is thus 

i \ 1/2 

A<f>p t corefg&vir<f> KTa , v \ ( 171* f \ ( M< 

c P e core < Ap 



1 /2 

,„ 110C = 2.50 x 10- ( ^#^ 1 (j^tr) EfM yr-, (54) 



3OM y V4000M Q 



_ 7.9 x 10- 4 ( ^] ( ^L.) ( Mcl ") S 3 / 4 M yr- 1 . (55) 
V 1.6 J \3OM ) \4000M Q ) c1 ° y y ' 

Thus we see that a massive core will tend to interact with clump material at a rate that is comparable to the rate at 
which it is collapsing or being eroded by star formation, m*/e COIC . This rate is a factor of several smaller in super star 
clusters with M c \ — 10 6 M o . 
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What is the fate of this swept-up mass? Gas between the cores ("intercore gas") is likely to have an AlfVen velocity 
that is greater than the Keplerian velocity at the surface of the core in view of its low density; it will therefore not be 
bound to the core after being swept up. It is likely that most of the mass is in cores, however, and insofar as the core mass 
distribution follows the IMF, most of this mass will be in low-mass cores. When one of these cores strikes the outer region 
of a massive core, it is quite likely that it will strip gas from the massive core, since the relative velocities are ^> v csc . On 
the other hand, if it strikes the inner parts of the massive core, it is likely to be absorbed by the larger core. Accretion of 
ambient material onto the massive core will decrease as soon as a protostellar outflow forms, and it will decrease further 
when the protostar becomes massive enough to generate significant radiation pressure. Thus, the interaction of the core 
with the ambient clump will result in an exchange of material that most likely causes the massive core to grow, but at a 
smaller rate than suggested by equation (55). In fact, this growth may be the mechanism that causes an initially bound, 
but stable, core to evolve into an unstable core that collapses to form a star. This growth will not be a runaway process 
since i? C ore cannot grow significantly once the expansion wave reaches the edge of the core. 

4.2.2. Bondi-Hoyle accretion by the star 

In the case of competitive accretion, it is the protostar that accretes gas from the clump. In this case, the appropriate 
accretion radius is given by i? acc = Rbh — 2Gm*/ (v 2 e] + c 2 ) (Bondi & Hoyle 1944). We generalize this result by including 
the effects of magnetic fields, c 2 = </>b(7 2 , so that 

n / <t>p, core/ <\>P, core \ 2Gto« 

^BH = ' - , 7 T3T- ( 56 ) 
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Evaluating this expression in terms of the virial parameter a v ir, we obtain 

10 / 4>p, corc/<M core \ / »^*_ 



where we made use of the relation R cl = (M c i/n'£ c i) 1 / 2 = 0.516(M cl /4000 M ) 1 / 2 E d 1/2 pc. 

This estimate for the accretion radius is actually an upper limit for two reasons: First, it is based on the assumption 
that the accreting gas is of uniform density (on scales of Rbh) so that it can efficiently dissipate in the wake of the star. 
Second, Bondi-Hoyle accretion will cease once the stellar mass is J> IOMq, since then the main sequence radiation pressure 
is enough to counter the star's gravitational attraction, thereby preventing the focusing of gas streamlines in the wake of 
the star (Bonncll et al. 1998). We therefore normalize our expressions to a stellar mass of 10 M in this subsection. 

With this result for the accretion radius, we find that the Bondi-Hoyle accretion rate in a star-forming clump is 



?ti*,bh = 1-66 x 10 



f 6 5/2 

J g'rp, core 



«vfr 2 4(Le(3 + 0B) 2 



10M Q J \ 4000M W J cl 



£ c rM yr-\ (59) 



- - 7 >< (iSj (ib^)" 5/4 ^ (60) 

We see that the Bondi-Hoyle accretion rate (for a 10M Q star, the maximum possible for BH accretion) is more than two 
orders of magnitude lower than the turbulent core accretion rate calculated in §4.1. If the clump survives for ~1 Myr while 
it forms a star cluster then the maximum fractional change in a star's mass is <~ 0.1(to*/10Mq)(M c i/4000M q )~ 5 / 4 £ c1 3/ ' 4 , 
so this process will not significantly alter the stellar mass under typical conditions of massive star formation. 

5. TWO-COMPONENT CORES 

Cores in molecular clouds are usually observed to have central regions dominated by thermal motions and envelopes 
dominated by nonthermal motions. The presence of a significant thermal pressure component in the inner regions of 
cores breaks the self-similarity of the structure observed on larger scales in molecular clouds (Williams et al. 2000) . As 
discussed in §3, the nonthermal pressure component is modeled with 7 P < 1 (Maloney 1988). Two components can be 
modeled as a multipressure polytrope in which P = P th + P n th, where P t h = P c th ^ s thermal pressure and P nt h is the 
nonthermal pressure (McKee & Holliman 1999). To allow for damping of nonthermal motions on small scales, the gas 
core can also be modeled as a composite polytrope, in which the central region is thermally supported and the envelope is 
supported by both thermal and nonthermal motions (Curry & McKee 2000). For analytic modeling, however, the simple 
"TNT" (thermal plus nonthermal) model introduced by Myers & Fuller (1992) and developed by Caselli & Myers (1995) 
is convenient. In this model, the density is assumed to be of the form 

H7f + (7) 2 - < 61 > 

with 1 < k p < 2. We assume that (-R C orcAo) 2 ~ fcp ^ 1) so that the first term dominates at the surface of the core and the 
model is equivalent to the sum of a singular polytropic sphere and a singular isothermal sphere: 



Rr 



14 




Fig. 2. — Evolution of m» (top panels) and m» (bottom panels). The straight lines show the purely nonthermal model, while the curved 
lines show the two component (thermal [100K] + nonthermal) core model described in §5. The left panels show only m,^ = 30Mq models. 
The different cases are k p = 1.5 with S c i = 0.316, 1.0, 3.16 g cm -2 (short- dashed, solid, long-dashed, respectively), and k p = 1.75 with 
S c i = 1 g cm -2 (dotted line). The right panels show cases with m*f = 7.5, 30, 12OM0 (bottom to top, respectively) for S c i = 1 g cm -2 and 
/.•„ = 1.5. 




Fig. 3. — Star formation time t t f versus final stellar mass m*/ for accretion from cores with k p = 1.5 at E G ] = 0.316, 1, 3.16 g cm 
{short- dashed, solid, long-dashed lines) and with k p = 1.75 at E c i = 1 g cm -2 (dotted line). As in Figure 2, the straight lines show the purely 
nonthermal model, while the curved lines show the two component (thermal [100K] + nonthermal) core model described in §5. 
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where, by construction, p s is very nearly equal to the density at the surface; it is therefore given by equation (8). 

To determine the mass-accretion rate m* and the star-formation time we need to know the value of the free-fall 
time, iff oc [/?(m)] _1 / 2 . Since the TNT model is just a rough approximation to the actual density distribution, we can use 
equation (3) to replace the TNT density distribution p(r) with a TNT density distribution p(M), 

2(i-j) 



P = Ps 



M 



+ 



M 



2(l-j) 



Ac th 
7rG 3 M 2 ' 

M th 
M 



2j 



(63) 



(64) 



where j is defined in equation (39). Here M tn is the mass below which the thermal density distribution dominates, and 
it will be evaluated numerically below. 3 
The protostellar accretion rate is 

\ 1/2 

(65) 

.Ps) 

from equation (28). Let 0* t h = 0.9757rv / 3/[8(l + Hq) 1 / 2 ] be the value of for the isothermal component, and 0» nt h be 
the value for the nonthermal component. Generalizing equation (65) with the aid of equation (64) for the density of a 
two-component core, we find 



to* = 



TO* TO* 
OC 
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TO*/ 



2j 
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th 



TO*/ 



2/ 



1/2 



(66) 



At early times, when the mass of the protostar is less than e C orcM t h, the accretion rate is constant — it is given by the 
isothermal rate (eq. 1) corrected by a factor e C orc/(l + Hq) 1 / 2 for the star formation efficiency and the effects of the 
magnetic field. Once the mass of the protostar exceeds this characteristic mass, the accretion rate rapidly approaches the 
nonthermal rate evaluated in §4. The evolution in time of the accretion rate and instantaneous stellar mass are shown in 
Figure 2 for two-component cores with a, T — 100 K thermal core. The star formation time versus final stellar mass is 
shown in Figure 3. 

We evaluate the mass of the thermal core, M tn , in terms of the Bonnor-Ebert mass, which is the maximum mass of a 
stable isothermal sphere of gas (Bonnor 1956; Ebert 1955), 
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where the numerical value is for a mean particle mass of 2.33toh- Equations (8) and (67) then imply 



M th = M BE 



32 



1.182 6 7T 3 ,4fcf, 
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(67) 
(68) 

(69) 
(70) 



If M coro is comparable to Mbe, then M t h — Mbe; as M corc increases, M t h decreases and becomes a negligible fraction of 
the core mass. Using equation (69), we can re-express the accretion rate as 



to* ~ TO*/ 



TO*/ 



3(2-fc p )/(3~fc p ) 



+ 



<t>* th 



0* nth 



32 



1.182 6 7T 3 ^fc| 3 



1/4 



BE 



TO*/ 



3/2' 



1/2 



(71) 



which makes the dependence on k p clearer. This form of the accretion rate shows that the thermal term dominates 



throughout the accretion only if to 



*/ 



< 



,M, 



BE- 



In the fiducial case, the condition for the thermal term to dominate 



throughout is to*/ < 1.5 x 10 (T/20 K) 2 /E c i M Q ; even for a temperature of 100 K, this is only a fraction of a solar 
mass. We conclude that the thermal core is not important in the formation of stars with to*/ > 1 M in typical regions 
of massive star formation. 



3 It must be kept in mind that M t h is defined in terms of the density profile and not the sound speed. Evaluating c 2 /c^ h with the aid of 
equation (6), we find 

4 \ 1/3 / M \ 2j/3 



For an isothermal gas (-y p = 1, k p = 2 and j = 0), we have c 2 = c 2 h as expected. However, as 7 P decreases below unity the value of c 2 at 
M = increases above c 2 h . Although the nonthermal pressure is greater than the thermal pressure at Mth, the two pressure gradients are 
comparable there. 
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Finally, we note that the fact that the pressure in regions of massive star formation is so high means that the minimum 
mass of an isolated star that can form is quite small, 

m* min ~ e core M B E- (72) 

If the temperature is of order 20 K, as in the case of clustered low-mass star formation (Jijina, Myers, & Adams 1999) 
and the pressure is that of a typical clump that is forming massive stars (P S;C oro ~ -Pel ~ 0.4 x 10 9 K cm~ 3 , corresponding 
to 4>p t C oro = 1), then this minimum mass is slightly less than 0.04M Q . As star formation proceeds in the clump, the 
temperature will rise and so will this typical minimum mass, m* m i n oc Mbe oc T 2 . On the other hand, pressure 
fluctuations in the clump can reduce the value of Mbe and allow the formation of stars with lower masses. The turnover 
in the IMF at low masses corresponds to the paucity of cores that form in regions of unusually high pressure and low 
temperature, where Mbe arid therefore TO* m ; n are small. In the Trapezium, the IMF has a broad maximum between 0.6 
M Q and the hydrogen burning limit (Muench et al. 2002). With E c ] = 0.24 g cm~ 2 from Tabic 1, (j>p yCO i C = 1 (since 
the low-mass stars are not concentrated in the center of the clump) and e coro = 0.5, equations (67) and (72) predict 
m» m i n ~ 0.15 M Q in Orion, which is qualitatively consistent with this observation. 

6. PROTOSTELLAR EVOLUTION, LUMINOSITIES AND COMPARISON TO HOT MOLECULAR CORES 

The properties of protostars depend on their accretion rates (e.g., Stahler, Shu & Taam 1980; Stahler 1988; Palla & 
Stahler 1992; Nakano et al. 2000). In particular, the accretion luminosity is 

= f acc Gm.m. = / / _m. rn* lO^A 

r* \0.5 J V 30M Q 10- 4 M Q yr" 1 r* J °' 1 ' 

where / acc < 1 is a factor that accounts for the energy advected into the star or used to drive protostellar outflows. If one 
were able to directly observe the protostar, it would be possible to distinguish the protostellar luminosity from that of the 
accretion disk. This is generally not possible for massive stars, so here we shall model the luminosities of hot molecular 
cores (HMCs) , 4 which include the contributions of both the protostar and the accretion disk. The advection of energy 
into the star at the accretion shock at the stellar surface is usually negligible. However, the mechanical luminosity of 
the protostellar outflow can reduce / acc below unity. In Paper II, we show that magneto-centrifugally driven protostellar 
outflows, which are powered by the rotational energy of the star-disk system, typically extract about half the total available 
luminosity from gravitational infall. Thus we set / acc = 0.5 and note that L acc ~ L w . 

For stars massive enough to be on the main sequence (to* > 20 M Q in regions of massive star formation — see below) , 
we can use the zero age main sequence (ZAMS) mass-radius relation to evaluate r„. From the solar metallicity models of 
Schallcr et al. (1992), we find r* = 7.27(m st /3OM ) - 55 R Q , provided to* > 1.5M Q . The accretion luminosity for stars on 
the main sequence is then 

p \ / \ 0.45 / 

/acc \ / TO* \ / TO* 



i ra = 6360(^— j ( JiRS _ T j L , (71) 

The final expression for L acc (eq.[75]) uses the fiducial massive star accretion rate from equation (41), again for protostars 
accreting on the ZAMS. For smaller masses the protostar has a size larger than the main sequence size, and the accretion 
luminosity is correspondingly reduced. 

In order to calculate the size of the protostar before it reaches the main sequence one must model protostellar evolution 
as the star accretes material and burns deuterium. D-burning provides a source of energy to support the star in a 
configuration with a central temperature ~ 10 6 K, so that the star is larger than when on the main sequence (central 
temperature of ~ 10 7 K). After a protostar has exhausted its internal store of D, the luminosity available from this process 
becomes limited by the rate of D accretion 5 : 



^2.5 x 10- 5 J V1O- 4 M yr" 

Once the star is so massive that this luminosity provides insufficient support, it starts to contract towards the main 
sequence. This takes of order a Kclvin-Hclmholtz time. If the star forms a radiative core during the contraction, shell 
burning of deuterium temporarily swells the protostellar radius by about a factor of two (Palla & Stahler 1992). 

Most previous investigations of protostellar evolution considered constant accretion rates (e.g. Stahler, Shu & Taam 
1980; Stahler 1988; Palla & Stahler 1992; Nakano et al. 2000), although Behrend & Maeder (2001) considered evolution 

4 In our terminology, some observed hot molecular "cores" should properly be termed hot molecular "clumps" . The examples chosen here are 
likely to be true cores. 

5 For consistency with the models of Nakano et al. (1995; 2000) and Palla & Stahler (1991; 1992) we have adopted an interstellar abundance 
of D relative to H of 2.5 X 10~ 5 (Bruston et al. 1981; Vidal-Madjar & Gry 1984). More recent observations tend to favor somewhat smaller 
abundances in the range D/H ~ 1.5 - 2.2 X 10~ 5 (Sonncborn et al. 2000; Moos et al. 2002; Vidal-Madjar & Ferlet 2002). These results also 
indicate that there are spatial variations in the D abundance within the Galaxy. Our adopted value is at the upper end of this observed range. 
Using a value of D/H = 1.6 X 10~ 5 reduces r„ by about 30% for protostellar masses > Mq (see also Stahler 1988). This raises Lacc and thus 
reduces the estimates of m* listed in Table 3. However, the effect is very small (< 10%) for m* ~ IOMq, becomes even smaller at higher 
masses as the internal luminosity becomes more dominant, and disappears completely once the protostar has joined the ZAMS. 



18 



with a growing accretion rate. To consider protostcllar evolution with an accretion rate evolving according to equation 
(66), we have developed a simple model based on that of Nakano et al. (1995; 2000). This model accounts for the total 
energy of the protostar as it accretes and dissociates matter and, if the central temperature T c > 10 6 K, burns deuterium. 
At all times the protostar is assumed to be a polytropic sphere of index, n, so that 7 P = 1 + 1/n. We have modified Nakano 
et al.'s model to include additional processes, such as deuterium shell burning and the transition between convective and 
radiative structure, and calibrated these modifications against the more detailed calculations of Stahler (1988) and Palla 
& Stahler (1991; 1992). 

The basic features of our protostellar evolution model, including its differences compared to the Nakano et al. model, 
are as follows. The initial conditions are set at to* = 0.1 M© and have a radius that increases with the initial accretion rate 
(Stahler 1988). Initially the star is contracting and supporting itself by the release of gravitational energy. At this stage 
convection is imperfect so that n is intermediate between 1.5 and 3 with a weak dependence on to*, which we calibrate 
from Palla & Stahler (1992). Once the star's central temperature reaches that required for D core burning, contraction 
halts and the star becomes fully convective (n = 1.5). We set the D core burning temperature at 1.5 x 10 6 K (Stahler 
1988). While there is sufficient D to burn, the star maintains this central temperature and its radius grows linearly with 
mass. After the protostar's D has been exhausted, D burning occurs only at the rate at which it is accreted into the star. 
Inevitably, as to* grows, the D luminosity (eq. [76]) becomes insufficient to support the star and contraction continues, 
raising the central temperature. 

A radiative barrier forms when the luminosity that can be transported radiatively, L ra( j, rises to equal L^. In the Nakano 
et al. model L ra( j is approximated by the main-sequence luminosity, L ms , of a star of the same mass. We generalize this 
and assume that a radiative barrier forms when Ld < /rad^ms, where for our adopted zero age main sequence (ZAMS) 
(Schaller et al. 1992) and by comparison to Palla & Stahler (1992), / rac j = 0.33. We approximate the effect of the 
formation of the radiative barrier by changing the polytropic index of the whole star from n — 1.5 to 3. Also at this 
point D can no longer penetrate to the stellar center so D core burning ceases. However, D will still continue to burn in 
a shell in regions hotter than ~ 10 6 K. We approximate the onset of shell burning to be simultaneous with the formation 
of a radiative barrier (see Palla & Stahler 1991 for a more precise description) and we assume it burns at the rate it is 
accreted. Shell burning swells the stellar radius by approximately a factor of two (Palla & Stahler 1992), which we take 
to be a constant in our models, equal to 2.1. In reality the amount of swelling will decrease as to* increases, but we 
find our simple approximation to be accurate to within <~ 10% compared to Palla & Stahler's results for this stage in the 
evolution. When T c £ 10 7 K, some hydrogen burning reactions begin and cause the star to become convective once more. 
The star is rapidly shrinking at this stage, and once its size becomes equal to the ZAMS radius (Schaller et al. 1992), we 
assume the protostar follows this mass-radius relation for the rest of its growth. 

We find our simple one-zone protostellar evolution model matches the mass-radius relations of the 10~ 5 , 3 x 10~ 5 and 
10~ 4 M© constant accretion rate cases with photospheric boundary conditions (no advection) of Palla & Stahler (1992) to 
better than about 10% over the entire evolution to the main sequence. Furthermore this range of accretion rates is similar 
to that predicted by equation (41) for values of to* relevant to pre-main-sequence evolution and for densities typical of 
high-mass star-forming clumps. The simplicity of our model allows us to apply it in numerical studies of the formation of 
clusters of many stars (Tan & McKee 2002b). 

Figure 4 shows the evolution of protostellar radius for 30 M© stars forming from cores with k p = 1.5 and 1.75 and 
at pressures equivalent to E c i = 1 gcm~ 2 . The accretion rate, here and for all cases in this paper, is derived using the 
two-component core accretion model (equation 66) with T = 100 K. As already discussed, the effect of the thermal core 
is quite small for the high-pressure conditions considered here. The accretion rate is dominated by gas supported by the 
nonthermal component. The smaller accretion rate of the k p = 1.5 case leads to ~ 25% smaller protostellar radii during 
most of the evolution compared to the k p = 1.75 case. Figure 5 shows the protostellar radii of stars of different final 
masses and illustrates the effect of increasing and decreasing the core's surface pressure by factors of 10. From equation 

(41), the accretion rate increases as n%f 4 , (to*/to*/) 5 (for k p = 1.5) and P s 3 ^ 8 oc E^ 4 . Accreting protostars join the 
main sequence at higher masses for higher accretion rates (Palla & Stahler 1992). For our fiducial case (k p — 1.5 and 
S c i = 1 g cm -2 ), stars join the main sequence at about 20 Mq. For larger values of k p or E c i, the accretion rates are 
higher and stars join the main sequence at somewhat larger masses (e.g., for to*/ = 30 M©, a star joins the main sequence 
at about 24 M© if either k p = 1.75 or E c i ~ 3 g cm -2 ). The maximum protostellar radius also increases with accretion 
rate and thus E of the clump. Again for k p — 1.5, 30 M© stars reach 10-14 i?© for values of E c i within a factor of 3 of 
those of the typical Plume et al. clumps. 

Our model for the evolution of the protostellar radius allows us to make predictions for the evolution of the bolomctric 
luminosity of embedded protostars (Figure 6). The bolometric luminosity is the sum of the internal and accretion 
luminosities, Lboi = L int + L acc . We have followed the approximation introduced by Nakano et al. to set the value of L- mt 
approximately equal to L ms , the main sequence luminosity corresponding to the particular value of to*, and we evaluate 
L ms from the ZAMS models of Schaller et al. (1992). Note that before the protostar joins the ZAMS, the source of L- lnt 
is gravitational contraction and D burning, rather than H burning. This approximation is least accurate at small stellar 
masses, but in this regime the total luminosity is dominated by the accretion luminosity. In fact, it should be noted that 
in the fiducial conditions of massive star forming regions, the turbulent core accretion model predicts quite high accretion 
luminosities (up to ~ 10 3 L©) when massive protostars are in their early (low-mass, <J 3M©) stages of formation. 

We compare our theoretical tracks with the observed bolometric luminosities of several HMCs thought to be illuminated 
internally by massive protostars. For the purpose of this simple analysis we consider three values for the initial column 
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density of the clumps in which these stars started forming, S c i = 0.32, 1.0 and 3.2 g cm -2 . This illustrates the effect of a 
change in the core's surface pressure by a factor of hundred. We also consider values for k p — 1.5 and 1.75. The observed 
luminosity and its uncertainty constrain the properties of the protostar. The constraint on to* is shown by the extent of 
the horizontal error bar for four sources in Figure 6 for the S c i = 1 gem -2 , k p = 1.5 case. The full list of constraints (on 
to* , to*/ and to*) and source properties is shown in Table 3. 

If the pressure (or S c i) is lowered, then the stellar mass range implied by a given luminosity moves to higher masses, 
since the accretion luminosity is reduced and so must be compensated for by an increase in internal stellar luminosity. 
Similarly, if the pressure is raised, then the stellar mass constraints are lowered. For Lboi <; few x 10 4 L Q and for the 
conditions considered here, the internal stellar luminosity (essentially the main-sequence luminosity) tends to dominate 
that of accretion. In these cases the strong dependence of L ms on to* means that large pressure uncertainties create only 
small uncertainties in to*. 

Osorio et al. (1999) have modeled the infrared and submillimeter spectra of accreting protostars and their surrounding 
envelopes, and applied their model to the same sources shown in Figure 6 and Table 3. We note that uncertainties in 
the structure of the gas envelope and the possible contribution from additional surrounding gas cores or diffuse gas will 
affect the observed spectrum. Comparing results, our inferred stellar masses are similar, except for the lowest luminosity 
source G34.24+0.13MM; this discrepancy is a result of low value Osorio et al. adopted for the main-sequence luminosity 
of a 10 Mq star. We systematically favor smaller accretion rates by factors ^2 — 5. Their high accretion rates of 
<~ 10~ 3 M Q yr _1 for stars with to* ~ 10 M Q are difficult to achieve unless the pressure is increased substantially; for 
example, if the stars are all destined to reach to*/ ~ 30 Mq, pressure increases of a factor ~ 40 are required. 

7. CONCLUSIONS 

We have developed a model (that we term the Turbulent Core Model) for the formation of massive stars, which is 
an extension of the classic paradigm of low-mass star formation (Shu, Adams, & Lizano 1987) and is to be contrasted 
with models involving competitive accretion (Bonncll et al. 1997, 2001) or stellar collisions (Bonnell et al. 1998). The 
principal motivations for the latter models are the short formation timescales (and correspondingly high accretion rates) 
mandated by observations of short star cluster formation times (Palla & Stahler 1999) and theoretical considerations 
of radiation pressure feedback (Wolfire & Cassinelli 1987). Such accretion rates are difficult to justify in the standard 
picture of isothermal core collapse (Shu 1977). Collapse from cores with nonthermal pressure support can involve faster 
accretion rates (Stahler et al. 1980), but there has been no self-consistent theory for predicting both the normalization 
of the expected accretion rates and their evolution. This has led to a vast range of massive star accretion rates (10~ 6 — 
10~ 2 A#q yr _1 ) being considered in the recent literature (e.g. Bernasconi & Maeder 1996; MP97; Nakano et al. 2000). 
Conventional theories of massive star formation face further problems: massive stars form preferentially in the centers of 
stellar clusters (Bonnell & Davies 1998) where the crowded environment makes it difficult to understand the existence of 
massive pre-stellar cores, and the high densities and pressures lead to a small thermal Jeans mass that is only a fraction 

Table 3 

Properties of Accreting High-Mass Protostars 



Source [Ref.] a Lboi k p S c i m* m*/ to* 

(1O 3 L ) (gem- 2 ) (M ) (M ) (lO^Meyr- 1 ) 



W3(H 2 0) [1] 


10- 


-90 


1.5 


1 


9.8- 


- 23.4 


> 10.5 


2.1 - 


- 5.6 








1.5 


0.316 


10.7- 


- 25.1 


> 10.9 


0.9- 


- 2.5 








1.5 


3.16 


6.2- 


- 21.1 


> 9.5 


4.6- 


12.7 








1.75 


1 


7.4- 


- 23.4 


> 10.7 


2.1- 


- 7.7 


Orion hot core [2,3J 


12.6 


- 50 


1.5 


1 


10.9- 


- 18.1 


> 11.3 


2.2- 


- 5.0 






1.5 


0.316 


11.5- 


- 19.5 


> 10.7 


1.0- 


- 2.2 








1.5 


3.16 


6.9- 


17.1 


> 10.4 


4.9- 


11.4 








1.75 


1 


10.5- 


- 18.5 


> 11.6 


2.2- 


- 7.2 


IRAS 23385+6053 [4J 


8- 


32 


1.5 


1 


8.9- 


- 15.3 


> 9.8 


2.0- 


-4.6 






1.5 


0.316 


10.0- 


- 16.3 


> 10.1 


0.9- 


- 2.0 








1.5 


3.16 


6.6- 


- 14.7 


> 8.5 


4.2- 


10.5 








1.75 


1 


6.6- 


- 15.8 


> 10.0 


2.0- 


- 6.9 


G34.24+0.13MM [5J 


1.6 - 


-6.3 


1.5 


1 


3.2- 


-9.0 


> 4.2 


1.1 - 


- 3.4 






1.5 


0.316 


4.3 - 


-9.5 


> 6.0 


0.6- 


- 1.5 








1.5 


3.16 


2.2- 


- 7.6 


> 3.2 


2.0- 


- 6.3 








1.75 


1 


2.6- 


- 9.2 


> 4.3 


1.1 - 


- 5.2 



a Rcfcrcncos: (1) Wyrowski ct al. (1999); (2) Kaufman et al. (1998); (3) Wright et al. (1992); (4) Molinari et al. 
(1998); (5) Hunter ct al. (1998). 
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Fig. 4. — Protostellar radius versus m* for stars with m t j = 30 Mq and for E c i = 1 gem 2 . The solid and dashed lines show the cases 
with k p = 1.5 and 1.75, respectively. The dotted line shows the zero age main sequence radius from Schaller et al. (1992). 
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Fig. 5. — Radii of accreting protostars with k p = 1.5. Solid lines show stars of final mass 7.5, 15, 30, 60 and 120 Mq (bottom to top) 
accreting from cores embedded in a E c j = 1 g cm~ 2 clump, typical of Galactic regions observed by Plume et al. (1997). The long dashed and 
short dashed lines show a 30 Mq star forming in a clump with mean pressure 10 and 0.1 times this value, respectively. The dotted line shows 
the zero age main sequence radius from Schaller et al. (1992). 
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Fig. 6. — Properties of accreting protostars. Solid lines show stars of final mass 7.5, 15, 30, 60 and 120 Mq accreting from cores with 
k p = 1.5 embedded in a S c i = lg cm -2 clump, typical of Galactic regions observed by Plume et al. (1997). The long dashed and short dashed 
lines show a 30 Mq star forming in a clump with mean pressure 10 and 0.1 times this value, respectively. The dotted line shows the zero age 
main sequence luminosity from Schallcr et al. (1992). The luminosity step occurring around 5 to 8 Mq, depending on the model, corresponds 
to the onset of deuterium shell burning, which swells the protostellar radius by a factor of about two and thus reduces the accretion luminosity 
by the same factor. Four observed HMCs are shown. The vertical error bar illustrates the uncertainty in their bolometric luminosities. The 
horizontal error bar shows the corresponding range of allowed values of m* for the £ c i = 1 g cm -2 , k p = 1.5 models (see Table 3). 
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of a solar mass (Bonnell et al. 1998). 

The above criticisms, taken together with observational hints that massive stars may form differently from low-mass stars 
(central concentration in star clusters, high degree of equal-mass binarity, complex morphology and extreme energetics 
of outflows), have motivated collisional and competitive accretion models. However, these are also not without their 
problems: as noted by Bonnell ct al. (1998), competitive (Bondi-Hoylc) accretion is suppressed for stellar masses above 
~ IOMq because of radiation pressure feedback, while stellar collisions require extreme ;> 10 8 pc -3 stellar densities, which 
have never been observed. 

The Turbulent Core Model for massive star formation overcomes the difficulties of the standard accretion scenarios by 
incorporating the effects of the supersonic turbulence and high pressures observed in massive star-forming regions (Plume 
et al. 1997). The high pressures mean that cores that become unstable are necessarily very dense and small, leading to 
high accretion rates and no over-crowding. The turbulent and nonthermal nature of cores gives them substructure and 
thus the protostar's accretion rate will exhibit fluctuations about the mean. Note that while dumpiness in the core is 
an attribute of our model for massive star formation, it is not a requirement, as in the collisional model of Stahler et 
al (2000). The turbulence and nonthermal support of the clump (protocluster) gas determines the IMF for stars above 
about a solar mass, whereas the Bonnor-Ebert mass is important in determining the IMF at lower masses. While our 
theory does not aim to predict the IMF, we note that the empirical core mass function is not too different from the stellar 
one (Testi & Sargent 1998; Motte et al. 2001), at least up to ~ 5M Q , the maximum mass probed by these observations. 
More massive cores, with relatively simple, centrally-concentrated morphologies have been observed, e.g. the Orion hot 
core (Wright et al. 1992) and W3(H 2 0) (Wyrowski et al. 1999). The stellar IMF would then be set by the core mass 
function, modulated by e coro (m»). In a model in which e corc is set by the feedback from bipolar outflows, it is found to 
have a relatively weak dependence on stellar mass (Matzner & McKee 2000; Paper II). 

A natural question is why gravitational fragmentation does not continue down to the thermal Jeans (sub-solar) mass 
scale inside a massive core, leading to the formation of a cluster of low-mass stars instead of a single massive star. In 
the self-similar model we have presented, clumps and cores exhibit density fluctuations on all scales down to the thermal 
Jeans mass. We have assumed that clumps and cores are relatively long-lived (over at least several dynamical times), 
which implies that most of the density fluctuations existing at a particular instant are gravitationally stable. How good is 
this assumption? The age spread of pre-main sequence stars in the Orion nebula cluster is of order 1 Myr, which is about 
an order of magnitude greater than the expected free-fall time of a typical Plume ct al. (1997) clump. The approximately 
spherical morphologies of many clumps observed by Shirley et al. (2002) suggests that they have existed for at least a 
dynamical timescale. If the onset of gravitational instability is a relatively rare phenomenon in this environment (i.e., < 
20% of the mass of a clump is undergoing gravitational collapse at any time), as we have assumed, then it is perhaps 
not surprising that gravitational fragmentation in a collapsing core is unlikely. Gravitational fragmentation is further 
suppressed by the tidal field of the embedded stars and protostars. 

The main results of this paper are the following: 

1. The surface densities in observed regions of massive star formation and in star clusters that contain, or did contain, 
massive stars, are typically within a factor 4 of S c i ~ 1 g cm -2 . The corresponding mean pressures are P/k ~ CE 2 ^ ~ 
10 8 — 10 9 K cm~ 3 , much greater than in the diffuse ISM or the typical location in a GMC. 

2. Cores have column densities similar to that of the clump in which they are embedded. The mean density of a core 
significantly exceeds that of its natal clump, and the radius is less than the the tidal radius. 

3. Cores that form massive stars are supersonically turbulent; there is no need for the gas to become subsonic in order 
for star formation to occur. 

4. The star-formation time is several times the mean free-fall time of the core out of which the star forms, but is about 
equal to that of the region in which the core is embedded. 

5. The time for a massive star to form in a typical region of of massive star formation is about 10 5 yr, with a weak 

( m i/ 4 ) dependence on stellar mass and a somewhat stronger dependence on the surface density of the clump in which it 

is forming (£ cl 3,/4 ). This timescale is short compared to estimated cluster formation times, but long compared to the ages 
of observed supernova remnants, which sometimes have been invoked as star formation triggers. 

6. The corresponding accretion rate, approaching 1O~ 3 M yr -1 for the most massive stars, is high enough to overcome 
the radiation pressure due to the luminosity of the star. 

7. For the typical case we consider, in which the cores out of which the stars form have a density structure p cx r~ 1,5 , 
the protostcllar accretion rate grows with time as cx t. These density structures are consistent with observed clumps 
and cores, while in an appendix we have shown that logatropic models are inconsistent. 

8. The rate at which a core accretes mass from the ambient clump is comparable to the rate at which it processes 
matter into a star. Once the star has formed, subsequent Bondi-Hoyle accretion is negligible, particularly for massive 
stars. 

9. Presenting a calculation of the evolution of the radius of a protostar, wc determine the protostellar accretion 
luminosity. When the (eventually) massive protostar is still less than a few solar masses, this luminosity can be several 
hundred to a thousand solar luminosities. Massive protostars join the main sequence at around 20M Q . 

10. Application to observations of the Orion hot core suggests a current protostellar mass of between about 11 and 
18 Mq and an accretion rate of a fewxlO~ 4 Af yr -1 . Similar properties are estimated for W3(H 2 0), the Turner- Welch 
object. 

The incorporation of feedback, including protostcllar outflows, ionization and radiation pressure, is the subject of Paper 
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II. In particular the question of when feedback prevents accretion is addressed. The implications of this work for star 
cluster formation will be examined in a future paper (see Tan & McKee 2002b for an initial discussion). 
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APPENDIX 

A. MEAN PRESSURE IN A CLOUD 



We wish to determine the mean pressure in a cloud. We assume that the cloud is ellipsoidal, with radius R normal to 
the axis of symmetry and radius Z along the axis. The total mass of the cloud is the sum of the gas mass and the stellar 
mass, M = M g + M*. Recalling that the pressure is related to the effective isothermal sound speed c by P = pc 2 , we have 



V 



PiV -**>-(*sh)&< (A1) 



where p is the volume average of the gas density and (c 2 ) = M g 1 J c 2 dM g is the mass average of c 2 . 
The virial mass of a spherical cloud is defined by (e.g., Myers & Goodman 1988) 

>„ . - *{o 2 )R 
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where a is the one-dimensional velocity dispersion. If the pressure in the cloud is entirely due to thermal and turbulent 
motions of the gas, then a = c; if magnetic fields contribute, then a differs from c, and we define <f>B = (c 2 )/(a 2 ) > 1. 
The virial parameter for a spherical cloud is defined as (Bertoldi & McKee 1992), 

a . = - M^L (A3) 

Note that a v [ r is defined in terms of the total mass M; the gas mass is a fraction f g of the total, M g = f g M. 

These definitions are readily extended to non-spherical clouds (Bertoldi & McKee 1992). For an individual cloud, R is 
replaced by i? bs = (-RMax-Rmin) 1//2 j where -R ma x and i? m i n are the maximum and minimum observed radii in the plane of 
the sky. For an ensemble of clouds, the value of i? bs averaged over orientation is denoted by R m , and the virial mass and 
virial parameter are defined in terms of R m . In particular, axisymmetric ellipsoidal clouds with a volume 4irR 2 Z/3 have 
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in the range 1 < goO m < 1-41 for 0.3 < Z/R < 3. For both prolate clouds with Z = 2R and oblate clouds with Z = Q.5R, 
we have ^ ge om = 1.13. Recent observations of star-forming clumps do not show strong evidence for asphericity (Shirley 
et al. 2002). For our numerical estimates we consider the spherical case with <f) gcom = 1. 

Combining these results, we find that the mean pressure is related to the surface density of the cloud (including any 
embedded stars), £ = M/itR"^, by 

p_ ^ 37r/ g B (/) gcom Q; v i r ^ 

Note that this is an identity, since it follows directly from the definitions of the various quantities involved. The physics 
enters primarily in the determination of the virial parameter a v -„ , which measures the effect of gravity, and in the magnetic 
parameter <j>B- Below, we estimate that a v ; r ~ 1 and 4>b — 2.8. We define the factor in parentheses in equation (A5) as 
<j>p. If we set / g = 2/3 as in the text and ^ ge0 m = 1, then we have 

P = $pGY? ~ 0.88GS 2 . (A6) 



A.l. The Virial Parameter a v i r 

Existing data suggest that ot v -„ ~ 1 for gravitationally bound molecular clouds and cores. Solomon et al. (1987) studied 
GMCs in the first Galactic quadrant and found that the virial masses for the clouds they studied were close to those 
inferred from gamma ray observations. Such observations determine the CO to H2 conversion factor, X = N(H2)/Wco- 
The most recent determinations of X from EGRET data give X = 1.8 x 10 20 cm~ 2 /K km s _1 in the local ISM (Dame, 
Hartmann, & Thaddcus 2001) and X = 1.56 x 10 20 cm~ 2 /K km s" 1 in the Galactic plane (Hunter et al. 1997). 
Determining X from virial masses, Solomon et al. found A v ; r = 2.2 x 10 20 /Mg' 235 cm~ 2 /K km s _1 , where M§ is the 
cloud mass in units of 10 6 M Q . [Their masses have been increased by a factor 1.31 to allow for a reduction in the distance 
to the Galactic Center from 10 to 8.5 kpc (a factor 0.85), for an increase in the mass due to the extrapolation to zero 
intensity recommended by Solomon and Rivolo (1989; a factor 1.4), and a difference in the definition of M v -„ (a factor 
1.1).] Note that A v i r /A = N v -„/N = M v i T /M = a v i r - For a typical cloud mass of 10 6 M© in the Galactic plane, this gives 
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a v i r = X v i r /X (Galactic plane) = 1.4. Note that the gamma-ray determination of X is necessarily biased toward nearby 
clouds, which have an average mass somewhat less than 10 6 M , which could raise a v i r somewhat. In any case, Solomon 
et al. estimated their cloud masses based on the assumption that a v i r = 1.1; if it is in fact larger than this, then their 
masses should be reduced by a factor l.l/a v i r . Dame et al. (1986) surveyed the brightest CO clouds in the first Galactic 
quadrant and determined cloud masses with X = 1.6 x 10 20 cm~ 2 /K km s . Reducing their masses by 0.85 2 and their 
radii by 0.85 to adjust to our adopted distance to the Galactic Center, we find that the mean virial parameter in their 
cloud sample is a v i r = 1.3, with a standard deviation of 0.76. Their results do not give a dependence of a v i r on cloud 
mass as found by Solomon et al. (1987), although it must be kept in mind that the Dame et al. sample is considerably 
smaller. 

Both bound and unbound structures exist within GMCs. Bertoldi & McKee (1992) found that most 13 CO clumps 
are gravitationally unbound; however, the most massive ones, which were also the ones that showed evidence for star 
formation, appeared bound. Tachihara, Mizuno, & Fukui (2000) have studied the C 18 cores in Ophiuchus and shown 
that they are generally bound. In the northern region, the value of a v ir calculated from the mean properties of the clouds 
(and including the thermal contribution of H and He to a 2 ) is 2.2 (excluding one core with M V - 1T = 12M), whereas in 
the main body of Ophiuchus it is 0.9. In Taurus, Onishi et al. (1996) found 40 C 18 cores with a v ; r = 0.9 as calculated 
from the mean properties (and with the same correction to a v i r ). The fraction of these cores that have embedded young 
stellar objects is 16%, 33%, and 48% in the northern region, main body of p Oph, and Taurus, respectively. The two 
regions with the most active star formation thus have a v ; r ~ 1. However, these data refer only to regions of low-mass 
star formation. There are no data from which one can directly determine the virial parameter for massive-star forming 
regions such as those studied by Plume et al. (1997). 

The value of a v - lr can also be estimated from theory. Let T = |M(cr 2 ) be the total kinetic energy (including thermal 

energy) and W = |aGM 2 / R be the gravitational energy; the parameter a is the ratio of the gravitational energy to that 
of a uniform sphere. The virial parameter is then a v ir = 2aT/\W\. Including the effects of the surface pressure P s and 
the magnetic fields, the virial theorem implies 

= 7^7^(1-1^. (AT) 



(l-P s /P)\ 

where M is the magnetic energy (Bertoldi & McKee 1992). The value of .M/|W| is quite uncertain, however. 

In the text, we model the cores as singular polytropic spheres (SPSs). If we adopt the same model for the clump in 
which the cores are embedded, we have 

(a 2 ) [5(c 2 )R\ a SPS 

where 



2 \ 6 - 5j p J \ 7 P 

is the virial parameter of an SPS (see MP96 or McKee & Holliman 1999). As discussed in §3, observed clumps have density 
distributions characterized by values of k p ranging from 1 to 2, with 1.5 being typical. The value of 7 P corresponding to 
k p = 1.5 is 2/3, which gives asps = 15/4. Below, we estimate (j) B — 2.8 so that a v ; r = 1.34, which is comparable to the 
observed values. 

A. 2. The Magnetic Parameter cf) B 
The effective sound speed c 2 = P/p includes the contributions of both the gas pressure and the magnetic pressure, 

9 9 B 2 SB 2 , 

c = CT + ^ + 777— > A1 ° 
8irp 24np 

where a is the 1-D velocity dispersion, which can be inferred from observation. We have distinguished between the 
pressure associated with the background field B (pressure = energy density) and that associated with a turbulent field SB 
[assumed to have a random orientation so that pressure =(1/3) energy density]. Taking the mass average of this relation 
gives 

<-> = ^>Hi + 5t?)' < A "> 

where E B = J (B 2 /8tt)cIV is the background field energy, E K = |/ g M(a 2 ) is the total kinetic energy of the gas, and Eg B 
is the energy of the turbulent field. The second term can be rewritten in terms of an average Alfven Mach number tua, 
where 

o 3(cr 2 ) E K ,„ „ N 

m ^ id = E~ B 

We now show that the last term in equation (All) is about 1/3. Fluctuations in the magnetic field are associated 
with motions perpendicular to the field, and, for small amplitudes and negligible ambipolar diffusion, these motions are 
in equipartition with the fluctuating field (Zweibel & McKee 1995): {l/2)pv 2 _ = SB 2 /8tt. There are also motions along 
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the field; assuming that the velocity field is approximately isotropic (v\ = 2a 2 ), we find (2/3)Ek — Esb- Thus, for 
isotropic turbulence the last term is 1/3. One can show that the same result is obtained if one considers the fluctuations 
to be MHD waves in a cold plasma and evaluates the wave pressure using the results of Dewar (1970). This result is in 
reasonably good agreement with the simulations of supersonic MHD turbulence by Stone, Ostriker, and Gammie (1998), 
who found E SB /2E K = (0.28,0.30) for j3 = 8nP th /B 2 = (0.01,0.1), respectively. We therefore adopt E SB /2E K = 0.3 for 
low (3, supersonic turbulence, which is the type that appears to be characteristic in molecular clouds. The effect of the 
magnetic field is then described by 

03 = P) =1+ 2^ + 2^ = L3+ 2^ (A13) 

from equation (All). This equation does not apply to high-/? plasmas, which have <f>s — 1- 

The magnitude of the magnetic field in massive-star forming clumps is unknown. The most complete survey of magnetic 
field measurements in star-forming regions is that of Crutcher (1999). While this survey primarily deals with regions of 
low-mass star formation, it does include regions whose density is comparable to, and even greater than, that in the Plume 
et al. clumps. For regions in which the magnetic field could be measured, the median value of the Alfven Mach number 
is 1.0. The magnetic field factor is then <pB = 1-3 + 1.5 = 2.8. 

B. THE FAILURE OF LOGATROPES 

Turbulent regions have velocity dispersions that increase with scale. Larson (1981) found that molecular clouds have 
line widths that increase with size at a rate intermediate between that of Kolmogorov turbulence ((ct 2 ) 1 / 2 oc r 1 / 3 ) and 
that of Burgers' turbulence, which is appropriate for a system of shocks ((cr 2 ) 1 / 2 oc r 1//2 ). Subsequent observations have 
shown that GMCs (Solomon ct al 1987; Heyer & Schlocrb 1997) and low-mass cores (e.g., Caselli & Myers 1995) have line 
widths that increase approximately as r 1 / 2 . Larson (1981) also found that the density tends to fall off as 1/r (i.e., k p = 1), 
so that in these regions there is a relation between the velocity dispersion and density, (cr 2 ) 1 / 2 oc p -1 / 2 . As discussed 
in the text, regions of high-mass star formation typically have k p ~ 1.5, so they usually do not satisfy this line-width 
density relation. Nonetheless, inferences drawn about the properties of turbulence in low-mass cores and GMCs would 
have implications for high-mass cores as well. 

Before continuing, we note that Larson's relations primarily refer to properties of an ensemble of clouds, rather than 
to the structure of individual clouds. The relation n oc 1/r then means that the column densities of different clouds are 
all about the same, to within an order of magnitude. (This can be understood in terms of the conditions needed for the 
stability of molecular clouds-see Elmegreen 1989; McKee 1999.) If this relation applied within individual clouds, and if 
the line-width size relation (ct 2 ) 1 / 2 oc r 1 / 2 also applied, then the turbulent pressure p(cr 2 ) would be constant and would not 
exert a force. It is for this reason that the density that can be supported by a self-similar cloud in hydrostatic equilibrium 
(eq. 4) vanishes for k p = 1. 

Lizano & Shu (1989) pointed out that the line-width density relation (ct 2 ) 1 / 2 oc p -1 / 2 follows if one assumes that (1) 
the turbulent pressure is proportional to the logarithm of the density, P tU rb oc In p, since then the signal speed {dP/dp) 1 / 2 
varies as p -1 / 2 , and (2) the velocity dispersion ct varies as the signal speed (dP/dp) 1 ^ 2 . Including the thermal pressure of 
an isothermal gas, they wrote the total pressure as 

P = pc 2 h + K Hp/Po), (Bl) 
where if is a numerical constant and po is a reference density. Gchman ct al. (1996) termed this a "logatropic" equation 
of state. Note that logatropes are not self-similar, and therefore avoid the "no-force" problem that afflicts self-similar 
polytropes with k p = 1. 

MP96 pointed out a fundamental problem with this formulation: Since the central velocity dispersion is observed to be 
thermal, it is necessary to set po — p c , which Gchman ct al. (1996) did. However, if one then sets K ^> 1, as both Lizano 
& Shu (1989) and Gehman et al. (1996) did, then the pressure becomes negative when the density drops only slightly 
from its central value. For example, if K = 10, the pressure is already negative for p = 0.9p c . 

To overcome this problem, MP96 proposed an alternate form of the logatropic equation of state, 

P = P c [l + A MP ln(p/ Pc )}, (B2) 

where Amp is a numerical constant that they estimate is about 0.18 in observed molecular clouds. (We have added a 
subscript "MP" to the parameter A defined by McLaughlin & Pudritz in order to distinguish it from the parameter 
defined in the text. We emphasize that our model in the main text is based on the polytropic results of McLaughlin & 
Pudritz, not the logatropic ones.) They present a thorough study of the structure and stability of spherical clouds with 
this equation of state. In particular they show that for ^4mp -C 1 and at large r, the density approaches the form 

P = Pc(n/r), (B3) 

where 



^{^t) • <B4) 

indeed, the density is quite close to this value for r > 1.5n. Gravitationally bound clouds have radii R » r\. MP96 
focussed on clouds near the brink of instability, which have R/r\ = exp(l/AMP — 1/4) — > 200, where the numerical 
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evaluation is for Amp = 0.18. The collapse of logatropic spheres was studied analytically by MP97 and numerically by 
Reid, Pudritz, & Wadsley (2002). 

Unfortunately, this form of the logatropic equation of state is also fundamentally flawed: It satisfies the observed line- 
width size relation only when the cloud is gravitationally unbound; furthermore, it leads to a "droop " in the line-width size 
relation for clouds that are bound. To see this, we evaluate the properties of the MP-logatropc in the case r » n, so that 
p oc 1/r: 
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where the virial parameter a v i r is discussed in §A.l. Since the effects of magnetic fields have not been taken into account 
here, this value of the virial parameter corresponds to a non in the notation of MP96. If magnetic fields are important, 
then MP96 assume that a non is still given by the rightmost expression in equation (B7), while the observed value (which 
they denote a mag ) is reduced. In terms of the parameter <pB introduced in Appendix A, we have a mag = a non /<f)B- 
Evaluating the exponent in the line-width size relation for MP-logatropes, we find 

rfln<0 . 10 (Bg) 
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This result immediately reveals the two flaws described above. First, these logatropes obey the observed line- width size 
relation (a 2 ) oc r only when the virial parameter a v - lr 3> 10/3, corresponding to an observed virial parameter ct mag 3> 1; 
i.e., the clouds are unbound. Second, when clouds are strongly bound (i.e., nearly critical), the line width actually decreases 
with size. This effect is clearly apparent in Figure 2 of MP96. In particular, for a critical cloud, we have a v i r = 35/18 
(eq. B7; MP96) so that the slope of the line- width size relation is -5/7. Such a "droop" in the line- width size relation 
has never been observed. For example, the observations of four GMCs by Heyer & Schloerb (1997) show that in each 
case the velocity dispersion rises approximately as r 1 / 2 over the entire range of their observations. (However, it should 
be noted that observations of high-mass star-forming clumps are not yet sufficiently detailed to rule out the possibility of 
a droop in the line-width size relation in these objects.) We conclude that, where tested by observation, this form of the 
logatropic equation of state also fails. 
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